NASA Contractor Report 2947 


Insensitive Control 
Technology Development 




C. A. Harvey and R. E. Pope 


CONTRACT NASI-14476 
FEBRUARY 1978 





NASA Contractor Report 2947 


TECH LIBRARY KAFB, NM 


ilWilHIII 


□OblSSl 


Insensitive Control 
Technology Development 


C. A. Harvey and R. E. Pope 

Honeywell Systems & Research Center 
Minneapolis, Minnesota 


Prepared for 

Langley Research Center 

under Contract NASl-14476 


IM/\S/\ 

National Aeronautics 
and Space Administration 

Seienlifie and Technical 
Information Office 


1978 




TABLE OF CONTENTS 


Section Page 

I INTRODUCTION AND SUMMARY 1 

Finite Dimensional Inverse Concept 2 

Maximum Difficulty Concept 2 

II SYMBOLS 4 

Upper Case 4 

Lower Case 5 

Greek Symbols 6 

Upper Case 6 

Lower Case 6 

III FINITE DIMENSIONAL INVERSE COMPENSATOR .... 7 

Mathematical Formulation 9 

The Illustrative Example 16 

The C-5A Example 27 

Evaluation 34 

Conclusions and Recommendations 60 

rV INSENSITIVE CONTROL SYSTEM DESIGN VIA 

AN INFORMATION MATRIX APPROACH 6 2 

Mathematical Formulation 62 

Use of the Information Matrix 65 

iii 



TABLE OF CONTENTS (continued) 


Section Page 

Design Method 70 

Problem Solution 72 

An Expression for the Cost Functional 72 

Gradient Expressions 79 

Computational Algorithms 81 

General Overview 81 

Modifications to the Bartels -Stewart Algorithm • ... 82 

Solution of the Nonsymmetric Case 84 

Solution of the Adjoint Case 85 

Computational Requirements for J, V J, and M ... 87 

xC 00 

Iterative Algorithm for Finding K’i- 88 

Successive Substitutions Scheme 89 

Conjugate Gradient Method 90 

Preliminary Evaluations 92 

C-5A Controller Design Evaluation 101 

Design Approach 101 

Design Evaluation 108 

Comparisons and Conclusions 118 

V CONCLUSIONS AND RECOMMENDATIONS 120 


IV 



TABLE OF CONTENTS (concluded) 


Section Page 

APPENDIX A . . 122 

APPENDIX B 130 

REFERENCES • • 134 


V 



Figure 

1 

2 

3 

4 

5 

6 

7 

8 
9 

10 

11 

12 

13 

14 

15 

16 


LIST OF FIGURES 


Page 


Compensated System with Model 13 

Compensated System without Model 13 

On-Line Identifier 14 

Responses to Unit Step (a^ =0, o-g = -0.5) 19 

Error Ratios for Step Responses =0, “2 " "0-3) • • • 20 

Error Ratios for Step Responses = 0. 1, o-g = -0. 1) 

with 6 = 1 and Sensor Noise 22 

On-Line Identifier Responses to Step Input, u = 1 . . • . 24 

On-Line Identifier Responses to u = sin(lOt) 26 

Variational Output and Incremental Outputs for 

Unit Step fie., M Variation 36 


Response to u^ = 0. 02 for = 1. 0, M^^ = 0.8 ... 38 

Response to Ug = 0. 02 for = 1.0, M.^ =0.8 ... 39 

Response to u^ = 0. 02 for = 1.0, M.^^ = 0.8 

with recycle every 1.25 second 40 

Response to u. Doublet for q^ = cup = 1.0» = 0.8 . . . 42 

1 it 

Response to u„ Doublet for q„ = tUp = 1.0, M = 0.8 . . . 43 

z t t 

Gust Response for q^ = u>£ = 1.0, M^ =0.8 45 

Gust Response with Recycle for q^ = = 1.0, M.^ =0.8 46 


vi 


III III 



Ill 


LIST OF FIGURES (continued) 

Figure Page 

17 Response to u. = 0. 02 for q„ = 1 . 25, (u~ = 0, 75, 

M^^ = 0.80 ! 47 

18 Response to u„ = 0,02 for q„ = 1.25, = 0. 75, 

= 0.80 : 48 

Wf 

19a Case 4R and Case 5R' Responses to u. = 0.02 for 

= u)£ = 1. 0, = 0.8 51 

19b Case 4R and Case 5R' Responses to u„ = 0,02 for 

qf = = 1.0, = 0.8 52 

20a Case 5R' and Case 4R Responses to u = 0. 02 for 

q^ = = 1.0, «)^ = 0. 75 53 

20b Case 5R' and Case 4R Responses to u„ = 0. 02 for 

q^ = = 1.0, uj^ = 0. 75 54 

2la Case 5R' and Case 4R Responses to u = 0. 02 for 

q^ = 1.25, u,^ = = 1.0 55 

2lb Case 5R' and Case 4R Responses to u„ = 0,02 for 

qf = 1.25, ujf = = 1.0 56 

22 Case 5R' and Case 4R Gust Responses for 

q^ = 0)f = 1. 0, = 0.8 57 

23 Case 5R' and Case 4R Gust Responses for 

qf = = 1.0, toj = 0. 75 58 

vii 



LIST OF FIGURES (concluded) 


Figure Page 

24 Case 5R' and Case 4R Gust Responses for 

= 1. 25, yj^ = =1,0 ■ 59 

25 Sensitivity versus p ^ 

26 Cost Ratio versus 97 

27 Aileron Control Activity versus 

28 Standard Deviation of Uncertain Parameters versus . 99 

29 Case 4R Maneuver Load Performance 110 

30 Case 4R Gust Load Performance (Bending Moment) ... 110 

31 Case 4R Gust Load Performance (Torsion Moment) • • • 111 

32 Case 4R Short Period Frequency Ill 

33 Case 4R Short Period Damping 112 

34 Overall Relative Score Comparison 114 

35 Normalized Performance /Range Comparison 115 

36 Normalized Spec Violation Comparison 116 


viii 



LIST OF TABLES 


Table 

1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

11 

12 


Error Ratios for Various Values of or2J P 

Open-Loop Performance Comparison 

Response Vector and Quadratic Weights 

Closed-Loop Performance Corriparison 

Steady State Relative Errors 

Gust Response Statistics 

Incremental versus Linearized Partials Comparison . . . 

Parameter Estimate Standard Deviations and Other 
Pertinent Information for Different Gains . 

Feedback Controller Performance --Case 4R 

Nominal = 1. 0, u>£ = 1. 0, =1.0) 

Gain Matrix for Run #3A 

Information Matrix Controller Performance 

Evaluation Model--Case 4R 

Ranking of Insensitive Controllers Including the 
Information Matrix Approach 


Page 

23 

28 

29 

30 

49 

50 
100 

103 

106 

107 

109 

117 


! 

I 


IX 




SECTION I 


INTRODUCTION AND SUMMARY 


The work described in this report represents a continuation of the effort 
initiated under NASA contract NAS 1-13680, Study of Synthesis Techniques 
for Insensitive Aircraft Control Systems, The by-products of that contract 
were two new advanced theoretical concepts for insensitive controller design 
that had been developed by contract consultants Professor William A. 

Porter of the University of Michigan and Professor David L. Kleinman of 
the University of Connecticut. The concept developed by Professor Porter 
has been designated the Finite Dimensional Inverse method whereas that 
developed by Professor Kleinman has been termed the Maximum Difficulty 
concept. At the conclusion of the initial effort, neither of the concepts had 
been developed to a point where the resultant insensitive controller designs 
could be evaluated on a realistic flight control example. The objective, of 
this contract effort, NAS 1-14476, Insensitive Control Technology Devel- 
opment, was to extend the theoretical base of the two concepts to workable 
insensitive controller synthesis techniques and to evaluate the resultant 
insensitive controller designs on a realistic flight control example. As in 
the first study, the C-5A longitudinal dynamics model was used as the test . 

ft 

bed for evaluation. 

The results of the study of the two concepts are summarized below 


FINITE DIMENSIONAL INVERSE CONCEPT 


The present formulation of this concept is much more suited to trajectory- 
type sensitivity problems than to stationary flight control problems. The 
controller designed for the C-5A example using this concept involves time- 
varying gains. The controller's performance for this example fulfilled the 
theoretical predictions, but the present formulation of the concept limited 
the reduction in sensitivity to certain selected outputs. For the C-5A 
example design, the control surface displacements were suitable as 
selected outputs. Unfortunately, reductions in sensitivity of these outputs 
did not yield reductions in sensitivity of other important responses, such 
as bending and torsion moments. 

Despite these deficiencies with respect to the C-5A example, the concept 
exhibits- certain promising aspects such as on-line parameter identification 
and sensitivity reduction of minimum phase input-output systems. The 
study also indicates that although the computational requirements are 
severe, they are not beyond current capabilities. 

MAXIMUM DIFFICULTY CONCEPT 

The objective of the Maximum Difficulty concept was extended to exploration 
of a technique devised by Professor Kleinman which utilizes the Informa- 
tion Matrix element of parameter identification theory. In its complete 
form, this concept requires insensitive controller design for a flight 
condition with minimal controllability index and, at the same time, desen- 
sitizing system responses to variations in uncertain parameters. The 
latter function, involving the Information Matrix approach, received a 
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majority of the attention in the contract and will be the main subject of 
the Maximum Difficulty concept discussion. 

In addition to the development of the Information Matrix concept. Professor 
Kleinman participated heavily in the actual design and preliminary evalu- 
ations. Formal evaluation of the concept was performed on the 15 -state 
Case 4R residualized C-5A model using the criteria defined in the initial 
phase. The evaluation revealed that the insensitive controller designed 
with the Information Matrix approach performed as well as the top-ranked 
controllers of the previous study. 
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SECTION II 


UPPER 

B 

D, D(or) 
F(a) 
G(or), G 
H(a) 

I 

J 

K 

K 

o 

K(t) 

M 

M 

w 

N(a,b) 

N(t) 

Q 


SYMBOLS 


CASE 

Bending moment 

Control input coefficient matrix in response equation 
Plant coefficient matrix 
Control input coefficient matrix 

Coefficient matrix of state vector in output or response 
equation 

Identity matrix 

Performance index 

Gain matrix 

Nominal gain matrix 

A Riccati matrix 

Finite dimensional inverse of T , Information Matrix 

o' 

Stability derivative (pitching moment due to vertical 
velocity) 

Normal distribution (a = mean, b = standard deviation) 
Inverse of K(t) 

Weighting matrix 
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s 


Weighting matrix 


Sm( p) Sensitivity index 

T Torsion moment 

T(or) Input-output transformation 

Nominal plant input-output transformation 
W Weighting matrix 

X(t) Coefficient matrix in the differential equation 

representation of M 

LOWER CASE 


e. 

1 

p(t) 

q 

q(t) 

u 

V 

w 

X 

^i 

z 


Basis vectors 
State vector for M 
Dynamic pressure 
Alternate state vector for M 
Control input vector 
Input to M 
Output vector 
State vector 

Input vectors to the nominal plant giving outputs 
Output of M 
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Ill I 


GREEK SYMBOLS 
Upper Case 

A(t) Scaling matrix in the differential equation representation 

of M 

TT(t) Symmetric matrix function of the 5^(t) 

I 

Lower Case 

a Vector of parameters 

P Scalar design parameter in insensitive compensator 

implementations 

Scaling factors 

fiT Perturbation transformation 

6a Aileron displacement 

6e^ Inboard elevator displacement 

Q Damping ratio 

T| Scalar white noise input 

Outputs corresponding to specific inputs and plant 
variations 

p Noise/signal ratio 

a Standard deviation 

tu Frequency 
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SECTION III 


FINITE DIMENSIONAL INVERSE COMPENSATOR 


One of the major objectives of this study was the development and quantita- 
tive evaluation of the finite dimensional inverse technique for compensation 
for parameter uncertainty. The finite dimensional inverse concept, 
conceived in the previous insensitive controller study, is based on the concept 
of the a priori construction of a set of inverse functions which are derived 
from a finite number of input -output pair relationships. The input -output 
pairs are specified by type of input, selected output, and a combination of 
uncertainties that represent variations in plant behavior from a nominal or 
no -uncertainty condition. 

In operation, the measured outputs of the plant are used to determine the 
degree of mapping of plant outputs on the prestored outputs at off-normal 
conditions. The degree of mapping then dictates the formation of the 
feedback signals using the inverse functions that are used to compensate for 
the plant operating at off-nominal conditions. In essence, the inverse 
functions represent the change in control which is necessary to negate the 
effect of parameter uncertainties. 

In investigating the finite dimensional inverse concept, the specific goal of 
this part of the study was to implement the concept for a simple 
illustrative example and for the C-5A example^ and to examine its 
performance and limitations via simulation. 
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The illustrative example is a first order system with two parameters. 

This example served two purposes. First, it provided a simple problem 
for purposes of debugging software and examining numerical aspects of 
implementation. Second, it permitted extensive analysis of the effects 
of nonlinearity with respect to the plant parameters; variation in the design 
parameter which governs the attainable degree of insensitivity; and sensor 
noise. 

The C-5A example provided a more realistic test for the concept. Addi- 
tional effects which could be examined with this example were those asso- 
ciated with authentic disturbances, unmodeled dynamics, authentic types 
of parameter uncertainty, and the choice of outputs of interest. 

Two forms of the compensator were examined. One included a simulation 
of the nominal as a model; the other excluded this model. In general, 
the performance of the compensators was consistent with theoretical pre- 
dictions: both configurations yielded reductions in sensitivity. The second 

configuration exhibited a tendency toward an initial high gain instability for 
a constant sensitivity gain. This could be alleviated with a time-varying 
sensitivity gain. The concept as implemented is based on the assumption 
that the system outputs are linearly dependent on the parameter variations. 
The effect of actual nonlinear dependence was found to be significant in the 
sense that the reduction in sensitivity for large parameter variations 
differed significantly from the linear theoretical predictions. But the 
sensitivity was reduced even for large parameter variations. 
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In the C-5A example it was found that the compensator provides reduced 
sensitivity in terms of the outputs of the system, but not in terms of the 
system states or other system responses of interest. To be of real bene- 
fit for such an example, this deficiency would need to be remedied. It 
was also found that the compensators did provide reduced sensitivity to 
gust disturbances and were not seriously affected by unmodeled dynamics. 

Details of the mathematical formulation and the experimental results are 
described below. 

MATHEMATICAL FORMULATION 

Consider a linear system represented by an input -output transformation 
T(g;) with O' denoting an r-dimensional vector of parameters. The system 
may be represented in state variable form as 

T(q'): X = F(or) X + G(a)u, x(0) = 0 (1) 

w = H(Qf) X 

where T(or) maps the input u into the output w. Let us assume without 
loss of generality that the nominal value of the parameter vector is zero 
and let T^ denote the nominal system. We assume that the dependence of 
T(cy ) on a is sufficiently smooth so that linearization about zero is an 
adequate model. If [e^, e^, . . ., is a basis for R , then the lineariza- 
tion of 6T(u) becomes 

6T(u) = S (2) 

i=l 
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where 


6T = T(a) - T (3) 

o 

Of = 2 O', e. (4) 

. , 1 1 
1=1 

= [T(ep -T^]u, i = 1,2, ...,r (5) 

We note that generally depends on u. 

For the moment, suppose that u is known a priori so that the may be 
computed. Let 

X = Linear Span 52»-**»5r^ ^ Range 6T(u) (6) 


For simplicity of discussion, we will assume a single input u. The exten- 
sion to more than one input is straightforward. Assuming that X is con- 
tained in the range of T^, there exist functions {y^, yg, . , . ,y^} such that 

q = Toyf. i = 1,2, ...,r (7) 

We call a map M a finite -dimensional inverse of if M is linear, bounded, 
and satisfies 

M|i = y., i=l, . , r (8) 

Assuming the set of functions {^^^(t)} are linearly independent on every 

finite interval [0, T], the map M may be realized in the following state 

2 3 

variable or differential equation form. ’ Let M map v into z. Then M 
is given by 
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z(t) = y(t) A(t) K(t) p(t) (9) 

p(t) = - |X(t) p(t) + A(t) |(t) v(t), p(0) = 0 (10) 

where 

y(t) = row (y^(t), ygCt), . . . , y^(t)) (11) 

|(t) = col (q(t), |^(t)) (12) 

A(t) = diag (X..), X.. = q(s)ds]'‘i (13) 

X(t) = diag(x„), x„ = |?(t)/J’J^ |^(s)ds (14) 

and K(t) is a symmetric r x r matrix. The matrix K(t) is the inverse of a 
symmetric matrix N(t) and satisfies the Riccati differential equation. 

K(t) = |{X(t) K(t) +K(t) X(t)} - K(t)TT(t)K(t) (15) 

where the elements of the matrices N(t) and lT(t) are defined as 

“ Jo (16) 

iry(t) = 5^(s)ds 5^(s)ds (17) 

The vector p(t) may be computed in the following manner as an alternative 
to Equation (10). Let us define the vector q(t) as the solution of the follow- 
ing differential equation: 

q(t) = 5(t)*v(t), q (0) = 0 (18) 

Then p(t) is given by the equation 

p(t) = A(t)q(t) (19) 
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since 


^ [A(t)q(t)] = A(t) q(t) + A(t) q(t) 

- -|X(t) [A(t)q(t)] + A(t) |(t)v(t) (20) 

The numerical solution for p(t) via Equations (18) and (19) was found to be 
much better behaved than the solution for p(t) via Equation (10) for t near 
zero. 

The finite dimensional inverse may be used as a compensator to reduce 
sensitivity. It may also be used in a parameter identifier mode. Two com- 
pensator configurations are shown schematically in Figures 1 and 2. The 
on-line identifier configuration is shown schematically in Figure 3, 

Notation in these figures will be used in the remaining discussion. In 
particular, we note that 

w is the output of the uncompensated system 

w^ is the output of the nominal system 

w is the output of the compensated system with model 

w is the output of the compensated system without the model 

Similar notation is used for the state vectors of these systems. For example, 
X is the state vector of the uncompensated system. 
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Figure 1. Compensated System with Model 



Figure 2. Compensated System without Model 






Figure 3, On-Line Identifier 

The map M for the compensated systems was computed according to the 
equations given above with the [5^} precomputed. The only difference in 
M between Figures 1 and 2 is in the inputs and outputs. In Figure 1, the 
input V is taken to be w - w^ and the output is called z. In Figure 2 the 
input to M is w and the output is called z. The map M for the on-line 
identifier was computed in a similar manner, except that y(t) is called y, 

A 

the r X r identity matrix. The were replaced by computed on line 
from the equations 

i = [— ] i = 1,2, ...,r (21) 

1 a . 

^ a = 0 

The output, z , is called a . 
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The sensitivity to variations in a of the system shown in Figure 1 may be 
represented by the index: 

lim _i 

SJP) = IHI-0 ||(T p) - T^MT(a) - T^) ‘|| (22) 

where T^(cy, p) is the closed -loop map from u to w. For a single-input/ 
single-output system, this index is the limit as |ja|| tends to zero of the 
ratio of the error between the compensated system and the nominal system 
outputs to the error between the uncompensated system and the nominal 
system outputs, i. e. , 

S (p) = lim ^ ~ (23) 

ll^lho w - w^ 

In the general case of an arbitrary M, 

Tm(a,p) = (I + p T(a)M)“^ T(c^) (I + PMT^) (24) 

and it follows that 

S^(P) = 1! (I + PT^M)'^|I (25) 

Thus, with M = T we would have 

’ o 

S^( p) = (1 + p) ^ (26) 

The M that was constructed above is an approximation to T^ ^ and in fact 
is equal to T^" on a finite dimensional subspace. Thus, we can expect 
that the compensator will exhibit reduced sensitivity. 
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For the system of Figure 2, defining a similar sensitivity index, 

S(p) = lim |)(T(a, 3 ) - T ) (T(a) - T )"^|| (27) 

where T(a, 3 ) is the closed-loop map from u to w. This leads to the same 
result as obtained for p), i,e,, 

S( 3 ) = ||(I+ 3 T^M)‘^|| 

Our major interest was to assess the utility of this concept to reduce 
sensitivity. 

As a by-product of our computational analysis, we also examined the utility 
of the concept for on-line parameter identification for the simple illustrative 

A 

example. The configuration is shown in Figure 3. In this case, M is a 

A 

finite dimensional inverse to the map T(u) from a to (T(q>) - T ) u = w - w . 

o o 

THE ILLUSTRATIVE EXAMPLE 

Several experiments were performed with the following scalar system 
assuming two uncertain parameters: 

T(Qf^, of 2 ): x(t) = -(1 + Qfj) X (t) + (1 + 0 - 2 ) u(t), x( 0 ) = 0 
w(t) = x(t) 

where the nominal values of q;^ and a 2 are zero. In addition to examining 
the performance of the two compensators (with and without an explicit model), 
the performance of an on-line identifier was analyzed. 
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The partial inverse for the compensators was constructed assuming for 
simplicity that the input u was a unit step. As such, the response of the 
nominal system denoted by (t) and the input-output pairs, [yj.(t), gj^(t)]. 


corresponding to perturbations in are: 

x^(t) = 1 - e"^ (28) 

y^(t) = -(1 - e"S (29) 

5 ^(t) = te'^-(l-e"S (30) 

y2<t) = 1 (31) 

?2(t) = 1 - e"‘ (32) 

The pairs [y^(t), |j^(t)] were computed from the equations 

5^ (0) = 0 (33) 

?2 = -?2 ^ “ 2 ’ ? 2 ^®^ ^ 

Yi * q - i * 1. 2 (35) 


with = ff 2 “ ^ ~ 1 • Equations (33) and (34) which define and 

the variational equations associated with the parameters a and Equa- 

tions (35) which define y and y_ are the nominal input-output relations, 

1 u 


The on-line identifier was constructed by choosing the output to be 


a 


Estimate of 
Estimate of or 2 


(36) 
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( 37 ) 


In this case the input -output pairs [y^, 5 ^^] are given by 

= [1 0], y2 = CO 1] 

” '5i -V h ' 'h 5i«» = 0 

with being conaputed on line for arbitrary inputs u(t). 

Experiments were conducted to examine the effects on the compensator 
performance of 

• Magnitude and direction of the parameter vector a, 

• Magnitude of the design parameter p, and 

• Sensor noise. 

Transient responses for a five -second interval were computed for various 
combinations of 01 and p with and without sensor noise. Noise -to-signal 
ratios less than or equal to one caused no significant changes in performance. 
Quantitative results for the variations in a and p with no sensor noise are 
summarized below. 

A typical response plot of the output w(t) is shown in Figure 4 for a step 
input. 

The responses of the compensated systems are closer to the nominal 
response than is the uncompensated response. The effect of the parameter 
p is also evident in Figure 4. The larger p yields less deviation from the 
nominal. The error ratios corresponding to the responses of Figure 4 are 
shown in Figure 5. Both these figures display a significant initial transient 
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OUTPUT 


r 



ERROR RATIOS 



Figure 5, Error Ratios for Step Responses (a = 0, a 9 ~ -0,5) 

1 ^ 
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for the compensator without the model. This transient is caused by the 
initial high gain on the system output. In the compensator with the model, 
this high gain is multiplied by the error rather than the output, and the 
initial transient is greatly subdued. A time -varying p could be introduced 
to alleviate the initial transient of the compensator without the model. 

The theoretically predicted values for the error ratio is (p + 1) ^ which 
gives 1/2 for p = 1 and 2/3 for p = 1/2, These values are somewhat less 
than the steady state values shown in Figure 5. The parameters for this 
figure are = 0 and ^ 2 ' -1/2. The error ratios for a case with = 0 
and Q-g * -1/10 are shown in Figure 6 for jS * 1. In this case of smaller 
parameter magnitudes, the steady state ratios very closely approximate 
the theoretically predicted value of 1/2. The direction of the vector 
(cKj^, (^ 2 ^ also influences the degree of compensation. The error ratios 
after five seconds for 12 different values of(of., a^) and three values of p 
are given in Table 1. The data generally confirm the theory for small 
values of and and the predicted trend in p. The major deviation 
occurs when the parameters are equal and are of the same sign. In this 
case the steady state values of the outputs for the compensated and uncom- 
pensated systems are equal. Thus, the denominator of the error ratio is 
approaching zero, and the fact that the ratios have the magnitudes shown 
is an indication that the compensators behave well even in this case. 

The on-line identifier was also evaluated by computing five -second transient 
responses. Responses to a unit step input with zero initial condition for 
several values of (a,, a.^) are shown in Figure 7. Four general charac- 
teristics are evident in this figure. First, if of^ is zero, the estimate 
of ^2 is exact to within the computer word -length accuracy. Second, for 
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TIME (SECONDS) 


Figure 6. Error Ratios for Step Responses - 0.0, = 

with 0=1 and Sensor Noise 
















Figure 7. On-Line Identifier Responses to 
Step Input, u = 1 



nonzero or^, the estimate of ^2 is more accurate than the estimate of 

Third, for nonzero 0 ?^^, the accuracy of the estimate is greatest in the 

initial phase of the response. Fourth, accuracy is better for small mag- 

2 2 

nitudes than for large magnitudes of + 0^2 when of is nonzero. 


The first characteristic is a consequence of the fact that the system output 
is linear in a 2 zero initial conditions. This fact also contributes to the 
estimates of being more accurate than the estimates of a. . The third 
characteristic is a consequence of the fact that as the system approaches 
steady state, the estimate is given by 


r ^ 2 A. ^ ~ 



h 

/\ A 

h^2 
- 2 

a = 1/2 

ss ' 


Lhh 

^2 

ss 

_?2 . 


which gives the linear relation 

(“2 ■ = ^^“^o^ss 


(39) 


(40) 


For this example, the only correct estimate satisfying Equation (40) occurs 
when = 0 and ^2 arbitrary. The fourth characteristic is a manifesta- 
tion of the nonlinearity with respect to of the system output. 


Responses of the on-line identifier to a sinusoidal input are shown in 
Figure 8. The same general characteristic's as for the step input occur. 

Additional responses were computed for other values of (or^^, 

cases with "sensor noise" added in the simulation. The noise caused no 

serious degradation in performance. 
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In summary, the compensator and the on-line identifier fulfilled the theoret- 
ical expectations. 

THE C-5A EXAMPLE 

A seventh order model of the C-5A was used as the design model. This 
model was generated from the 79th order model described in Reference 1 
as Case 1. The unsteady aerodynamics were tmncated as in Reference 1 
to arrive at a 42nd order model called Case 2. Then a 13th order model 
called Case 5 was computed by the process of residualization. This model 
retained only one flexure mode. The seventh order gust model was then 
approximated by a first order gust model yielding the seventh order reduced 
model called Case 5R'. The data for this model are given in the Appendix A. 
The parameter variations considered were the same as those in Reference 1. 

In Reference 1, a 15th order model called Case 4R which retained two 
flexure modes was used as a design model. Comparison of the open-loop 
performance of Case 4R, Case 5R, and Case 5R' models is given in Table 2, 
The results indicate a high degree of consistency. 

A nominal controller for Case 5R' was computed using quadratic optional 
theory with the same response vector and weights as used in the nominal 
controller design for Case 4R. The response vector and quadratic weights 
are given in Table 3. The closed-loop performance for the nominal param- 
eter setting for Case 4R and Case 5R' is shown in Table 4. Again, there is 
a high degree of consistency. 
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TABLE 2. OPEN -LOOP 


Case 4R 

6 

Maneuver load 0. 0427 x 10 

bending (N-m) 

0 

Gust load bending 0* x 10 

(N-m) 

Gust load torsion 0, 179 x 10^ 

(N-m) 

Short period frequency 1.62 

(rad/ sec) 

Short period damping 0,56 

(sec ) 


3RFORMANCE COMPARISON 


Case 5R Case 5R' 


0.0427 X 10® 

0.0427 X 10® 

6 

0, 12 X 10 

0. 12 X 10® 

0.179x10® 

0. 174x 10® 

1,55 

1.55 

0.57 

0,57 


TABLE 3. RESPONSE VECTOR AND QUADRATIC WEIGHTS 


Response Vector 

Physical Quantity 

Weight 



= bending moment at wing root 

1 X lO”^^ 



= torsion moment at wing root 

1 X 10“^ 

% 


= rate of change of bending 
moment at wing root 

5 . 5 X 10 



= rate of change of torsion 
moment at wing root 

1 X lO"^^ 

-^5 


= aileron displacement 

0,32 X 10^ 


«e. 

1 

= inboard elevator displacement 

0 


6 

a 

m 

= function of aileron displace- 
ment and aileron command 

1 X 10® 


6 

e. 

1 

= inboard elevator rate 

1 X 10^ 

% 

^CF 

a control follower response 

2 X 10® 
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TABLE 4. CLOSED-LOC 


Maneuver load control 
bending % change 

Gust load — bending 
alleviation 

Percent change --torsion 
Handling u;^p(rad/sec) 

Qualities ^(sec 

Stability Gain 6a 
6e 


Phase 6a 
6e 


PERFORMANCE COMPARISON 


Criterion 

Case 4R 

Case 5R' 

< -30% 

-40% 

-41% 

< -30% 

-35% 

-34% 

< + 5% 

-31% 

-27% 

> 1.6 

2, 13 

2, 16 

0.7 - 0..8 

72 

0. 73 

> 6db 

00 

00 


29 db 

00 

> 45 

00 

00 


00 

00 
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The equations used to implement the compensator for this example are 
given below. In this case the input is two-dimensional and the inverse is 
based on independent step inputs in each channel. The equations are given 
for a general two-dimensional output. 

The equations for the nominal states and outputs are 


(F + G,K ) X. + G, u. 
o 1 O 1 1 1 

(41) 

(H + D K ) X. + D u. 
O O O 1 O 1 

(42) 


where i = 1, 2 represents step inputs on the aileron and inboard elevators, 
respectively, and the subscript o indicates nominal value. 

The equations for the variations in states and outputs are 

= (F^ +G^K^) X (i,j) +(F^ - F^) X. (43) 

I (i, j) = (H^ + D^K^) X (i, j) + (H. - H^) X. + (D. - D^) u. (44) 

where F., H., and D. are computed from the general expressions for F, 

3 3 3 

H, and D in terms of q„, tu~, and M with j = 1,2,3 denoting the following 

I I W£ 

specific variations: 

j = 1 : q^ = 1.0, = 1.0, = 0.8 

j = 2 : q^ = 1 . 0, (jDj = 0 . 75 , =1.0 

3 = 3 : q^ = 1. 25, (pj = 1-0, = 1.0 

The equations relating the variational outputs that are to be outputs of the 
nominal system are 
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( 45 ) 


X (i, j) = (F^ + G^K^) X (i, j) + y (ij) 
i(ij) =(H^+D^K^)x(i,j) (46) 

A 

where y(i,j) is the inverse input defined such that |(i,j) = |(i. j)* 
expression for y(i, j) will now be derived. This derivation will assume that 
the number of outputs equals the number of controls. 

Let us introduce 

e(i, j) = X (i, j) - X (i, j) (47) 

and 

e (ij) = 5(i,j) - i(i,j) (48) 

Then 

e (i, j) = Fq X (i, j) + (F^. - F^) x. - F^ x (iJ) - G^^y (i, j) 

= F^ e (i, j) + (F^. - F^) x^ - G^ y (i, j) (49) 

where 

F = (F + G, K ) 
o o 1 o 

Defining H = (H + DK ), Equation (48) may be written as 
o o o 

e (i, j) = e (i, j) + (H. - H^) x. + (D^ - D^) u. (50) 

Assuming u^^ = 0 and differentiating Equation (50), we have 
?(i,j) = H^e(i,j)+(H. -H^) X. 

= H^ [F^e(i,j)+(F. -F^)x. - G^ y (i, j)] + (H. - H^) (51) 
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Now y(i, j) is to be such that e(i, j) is identically zero, which is equivalent 

• A 

toe(i, j) being zero since x(i, j), x(i, j), x^^, and u^^ are initially zero. 
Thus, setting the left-hand side of Equation (51) to zero, we may solve for 
y(i, j) obtaining 

y(i,j) = C[H^F^e(i,j) + H^(F. -F^)x.] +(H. (52) 

Using this expression for y(i, j). Equation (49) may be rewritten as 
e(i,j) = [I - CF<,ea.i) + (P. - F^) x.] 

- Gj(H^Gj)'^ (H - H^) X. (53) 

Now the dot products (for i, j = 1,2,3) are defined as 


♦ 

V. . 

13 

= |(l,i) 

• ?(1.3) 

(54) 

ft 

V3,j+3 

= ?(2,i) 

• 1(2,3) 

(55) 

\j+3 


= ?(l,i) • |(2,j) 

(56) 

« 

‘’j 

= ?(1,3) 

• W J 

measured 

(57) 

* 

= 5(2,j) 

^measu red 

(58) 


Then defining 


“Ik “ ''kk “ 

, k = 1,2, . . .,6 (59) 


yU.j) if V.. f 0 

, j = 1.2,3 (60) 

0 if V.. = 0 
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^j+3 ' 

, j “ 1, 2, 3 

° ''j+3,j+3 " ° 

(61) 

K = 

N-1 

(62) 

11 

1, k = 1,2, . . .,6 

(63) 

I! 


(64) 


the output of the compensator is 

z = [71,72. •• -»ye] K p (65) 

A -1 

where y is yA and a = diag (• . 

In the first attempt to implement the compensator, the first two responses, 
bending moment and torsion moment, were chosen as the outputs. Unfor- 
tunately, the transfer matrix from the inputs to these outputs has a zero in 
the right half -plane. This caused the coefficient matrix in Equation (53), 
namely [I - Gj^(H^G^) ^ H^]F^, to have eigenvalues in the right half-plane 
leading to divergent functions y(i, j). Rather than implement this compen- 
sator with internal instability, it was decided to choose the aerodynamic 
surface positions, 6a and ee^^, as the outputs. This choice eliminated 
the internal instability. 

Evaluation 

Experiments were conducted to examine the quantitative performance of 
the compensator. The parameter variations considered were the three 
major parameter variations of Reference 1, i.e., dynamic pressure, q; 
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structural frequency, tu; and the stability derivative, Dynamic 

pressure variations cause changes in all the aerodynamic terms and, hence, 
most of the elements of the coefficient matrices. Structural frequency 
variations induce variations in a significant subset of coefficients in the 
complete model. The variation in permits examination of the effect 
of variation in a single coefficient. The physical motivation for treating 
these variations is discussed in Reference 1. This example also has per- 
mitted examination of the effect of authentic gust disturbances, A single 
sample on a five-second interval was used to examine the gust effect. 

To examine the effect of unmodeled dynamics, the compensator designed 
for the 7th order model was used in conjunction with the 15th order 
model. Experiments were also conducted to test the effect of the design 
parameter, p. The model was incoiporated in the compensator for all the 
C-5A experiments. Another experiment conducted tested the effect of re- 
cycling the time -varying gains of the compensator. These gains were 
computed off line for a five -second interval. 

In the recycle experiment, the gains for the first fourth of this interval were 
used repeatedly in each succeeding fourth of the interval. This experiment 
was motivated by two considerations. The first is that the example is 
essentially stationary and that infinite data lengths are impractical. The 
second is that the assumption of linearity of the outputs with respect to 
parameter variations is less valid for longer time intervals than for short 
ones. This phenomenon is shown in Figure 9. The two components of 
f(2, 1), the output of the variational equation associated with for an 
inboard elevator step input, are shown along with the actual increments 
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in the outputs corresponding to 4 and 20 percent variations in for the 
same step input. The nonlinear effect is clearly more pronounced in the 
later part of the five -second interval. This phenomenon is common to all 
the components of the |'s. 

To determine the effects of the design parameter, p, responses to gusts 
and step inputs of 0.02 radian magnitude in u^^ and Ug were computed. To 
determine the effects of different parameter variations, such responses 
were computed for independent variations of individual parameters and two 
cases of variations in all three parameters. These latter two cases were 
found to be "worst case" variations in Reference 1. 

Figures 10 and 11 show the deviations from nominal of control surface 
deflections in response to step inputs for p = 0, 0.5, 1.0, and 3.0 with the 
model in the loop and a 20 percent variation in The reduction in 

these deviations for increasing p is consistent with the theoretical pre- 
diction. 

The effect of recycling is shown with p = 1 in Figure 12 for the same param- 
eter variation and a step in u^. The recycling induces a severe transient 
following the start of each recycle. This is due in part to the fact that the 
compensator output is zero for an initial interval, and in part to the effec- 
tive high gain in the system following this zero output interwal. Although 
the recycling tended to deteriorate the 6a response, it improved the 
response. 

Figures 13 and 14 show the effect of changing the input. The inputs for 
these responses were chosen to be a positive step for 1 sec, followed by 





















a negative step for 1 sec, and then zero for 3 sec. This type of input is 
more realistic for an aircraft. 

The gust response for this parameter variation is shown in Figures 15 and 
16. The effect of p is shown in Figure 15, and the effect of recycling is 
shown in Figure 16, In Figure 15, the nonlinear effect is evident with 
better performance in the earlier portion of the interval than in the later 
portion. In this case, recycling seems to improve performance generally 
in spite of the induced transients. 

Figures 17 and 18 display the effect of p for step inputs for the so-called 
"worst" Case 1 condition. Again, the trend is consistent with the theory. 

The steady state relative errors for the parameter variations considered 
are summarized in Table 5 for step inputs , The trend is generally in 
accordance with the theory. Exceptions do occur in cases where one of 
the ratios is negative. 

Table 6 presents the gust response statistics for the same parameter 
variations. Here the mean and standard deviations are computed for the 
time series of numerical integration on the five -second interval for a 
single gust sample. Again, the results are generally in accord with the 
theory. 

The final experiment consisted of testing the effects of unmodeled dynamics. 
Figures 19, 20, and 21 present comparisons of step responses for the 
7th and l5th order systems with the three individual parameter variations. 
Gust responses for these cases are shown in Figures 22, 23, and 24. The 
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(6a - «a.)10 












Plant 


FI 




Step Input 


x(5) - X (5) 
0 



x(6) 

- X (6)’ 

0 

x(5) - X (5) 
0 J 

t = 5 


x(6) 

- X (6) 

0 J 


1 

3 

1 

3 


XL^ = 0.02 
= 0.02 
Ug = 0.02 
U 2 = 0. 02 


0.52 

0.27 

0.52 

0.27 


0.52 

0.26 

0.52 

0.26 


t = 5 


F2 1 
3 
1 
3 


= 0.02 
u ^ = 0.02 

U 2 = 0. 02 
U 2 = 0.02 


0.20 

-0.99 

0.48 

0.21 


0. 19 
0.30 
0.50 
0.23 


F3 1 
3 
1 
3 


= 0,02 
= 0.02 
U2 = 0. 02 
U2 = 0. 02 


0.56 
0.56 
0.49 
0. 25 


0.28 

-0.03 

0,65 

0.34 


WCl 1 
3 
1 
3 


u^ = 0. 02 
= 0, 02 
U2 = 0.02 
U2 = 0.02 


0, 27 
0. 70 
0.50 
0.24 


0. 22 

- 0.01 

0.87 

-0.04 


WC2 



0.02 

0.02 

0.02 

0.02 


- 0.01 

- 0.01 

0.57 

0.26 


0.34 
0.34 
0. 78 
0.71 
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Figure 23. Case 5R' and Case 4R Gust Responses for 







unmodeled dynamics have very little effect. The greatest effect appears to 
be the change in the gust sample caused by the added filtering. 

CONCLUSIONS AND RECOMMENDATIONS 

The performance of the compensator generally lives up to theoretical pre- 
dictions. There are four major areas where the current formulation is 
deficient for an aircraft flight control application such as the C-5A. They 
are 

1. Reduction of sensitivity of arbitrary responses of interest, 

2. Application to a stationary system, 

3. Adequate treatment of nonlinear dependence on parameters, 
and 

4. Severity of computational requirements. 

In the C-5A example, insensitivity of many responses and particularly 
bending and torsion moments are desired. Increases in these responses 
were caused by the compensator in most instances. Presumably, this 
would have been avoided if they were used as the outputs. But in this 
example, it would have required an internally unstable compensator. The 
formulation should be modified to include insensitivity to such responses. 

The recycling was an ad hoc attempt to modify the compensator to account 
for stationarity. It was not completely satisfactory; alternative formula- 
tions should be considered. 
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The nonlinear dependence on parameters is significant for the C-5A example 
and the theory should be modified to encompass this phenomena. In this 
example, this effect and the stationarity effect were coupled by the choice 
of inputs used. Other inputs or additional inputs could alleviate this coupling. 

The computational requirements associated with the implementation were 
significant for the seventh order design model with the limited number of 
inputs used in construction of the finite dimensional inverse. This was 
not a major concern in performing an evaluation of the concept, but it would 
be for an operational system. 

Thus, although the concept lived up to expectations, further development 
is required before the technique could lead to operational systems for an 
application such as the C-5A, 
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SECTION IV 


INSENSITIVE CONTROL SYSTEM DESIGN VIA 
AN INFORMATION MATRIX APPROACH 


The technique described in this section was developed by contract consul- 
tant Professor David L. Kleinman of the University of Connecticut, The 

4 

technique is based on the utilization of the Fisher Information Matrix 
which is a fundamental feature of maximum likelihood parameter identi- 
fication. In identification applications, it is desirable to minimize in some 
sense the inverse of the Information Matrix, or the dispersion matrix, in 
order to enhance the identifiability of a set of system parameters. With 
respect to sensitivity, given a set of responses, it was hypothesized that 
minimizing the Information Matrix itself would reduce the identifiability 
of system parameters and, consequently, the sensitivity of the system 
response to variations in those system parameters. 

The evaluation of control systems designed with the technique consisted of 
a preliminary evaluation of the effect of adjusting design parameters on a 
single system response of the C-5A test example and a full system 
evaluation for the 15th order model as was done in Reference 1. 

MATHEMATICAL FORMULATION 

Consistent with the assumptions given in Reference 1, the system to be 
controlled may be represented by a set of linear constant coefficient 
differential equations. 
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( 66 ) 


X = F(ff) X + G^u + Ggl] 


where 

X is an n state vector 

— X 

u is an n control vector 
— u 

T| is a scalar white noise with N(0, 1) 
a is an n^ parameter vector 

The n^ system responses may be represented by 

r = H(or)x + D(a)u (67) 

In Equations (66) and (67), the matrices F, G^^, Gg, H, D have the appro- 
priate dimensions. As in Reference 1, only F, H, D are assumed to 
depend on the parameters, a. However, G^ and Gg can also depend upon 
a, in general. 

The parameter vector is assumed to be bounded: 

Since the parameters are regarded as uncertciinty factors (as opposed to 
absolute deviations), the matrices F, H, and D about any operating point 
a can be expressed as 

p 

F(or) = F + E (a. - 1)F. (68a) 

— O . , 1 X 

1=1 

H(or) = H + S (a. - l)H. (68b) 

° i=l " 
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(68c) 


D(a ) = + 2 (or. - 1)D. 

i=l 

where F . H , and D are matrices at a nominal point, i. e., by definition 
o o o X- » 

where = 1. The matrices F^, etc., are the gradient matrices evaluated 
at the nominal point, 

F, 4 7,F(a) I . = jSl =grad (F) (69) 

= grad^ (H); = grad^ (D) 

For the system described by Equations (66) through (68), the problem is 
to determine a feedback control law u = Kx such that the system is "insen- 
sitive" to parameter variations and satisfies representative performance 
criteria. For this application, latter criteria are expressed via the 
minimization of a quadratic criterion 

= E{r^ Qr} (70) 

where the weighting matrix has been selected to meet specifications at a 
nominal design point. The minimization of yields a nominal feedback 
control 

u = K X (71) 

— o — 

Of course, this set of gains has been determined with the neglect of 
explicit sensitivity criteria. 
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Use of the Information Matrix 


The Fisher Information Matrix for the parameter set a and the responses 
r(t) in Equation (67) is given (approximately)* by 

“ = ® (i) 

where T is the observation time or measurement interval, S is a weighting 
matrix to be discussed later, and ftr/gjo is an n^ x p sensitivity matrix 
with elements 

\baj bof- 

ij 

Thus, the i-th column of the sensitivity matrix reflects the sensitivity of 
r with respect to the i-th parameter. 

The use of the information matrix for optimal input design to enhance the 
identification of unknown parameters from the measurement set r (assuming 
that an efficient, unbiased estimator of a exists) is well known. 

This technology has served as the motivation behind the present work. 

The major points relevant to the discussion are the following: 

A 

*The approximation used is to replace x by x, i.e, , we omit the Kalman 
filter innovations representation in the Maximum Likelihood formulation. 
This greatly simplifies the problem while retaining the essential param- 
eter sensitivities. 
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1. The inverse of the information matrix, called the dispersion 
matrix, gives a bound on the parameter error covariance 

E{(or - a)(a -a)'} ^ M ^ = (Cramer-Rao lower bound)^ 

Here, a is the actual parameter vector and a is the estimate. 
The idea behind optimal input design is to find u to maximize 
a metric on the information matrix (e.g., det (M), tr (M), 
max \ (M)) or equivalently to minimize a measure of the 
dispersion matrix. The net effect will be that one can place 
more confidence in the parameters estimated from the input- 
output data. 

2. Obviously, maximizing tr(M) would lead to increased 
sensitivity of r(t) with respect to the unknown parameters, 
inasmuch as M is a function of the output response sensi- 
tivities . 

With respect to sensitivity, the problem is precisely the 
opposite to optimal input design, i.e., to determine a control 
input such that r(t) is least sensitive to parameter variations. 
To solve this "inverse” problem, we seek the "worst" input 
from an identification viewpoint, i.e., one that makes the 
parameter set as unidentifiable as possible. Therefore, it 
is natural to seek to minimize a measure of the information 
matrix. 

Determining which measure to use poses a problem. In the 
identification literature, it is reported that choosing optimal 
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inputs to minimize the trace of the dispersion matrix gives 
the most accurate estimates of a. Therefore, it is conjec- 
tured that minimizing the trace, or weighted trace, of the 
information matrix should be the "best" criterion for choosing 
a desensitive feedback control. 


As the observation interval T -» <» in Equation (72), the norm 
of the information matrix approaches infinity. Thus, since 
we are dealing with a steady-state optimization problem, it 
is more appropriate to consider the information matrix per 
unit time, or the average information matrix, M^, In the 
steady state, we have approximately 



(73) 


The weighting matrix, S, is the inverse of the measurement 
noise covariance. Since the problem formulation does not 
include measurement noise (note that it could), S will be inter- 
preted in terms of a "pseudo" measurement noise injected 
onto r(t). We select S to be diagonal with elements 

Si=(pa^r^'^ i = l,...,nr (74) 

i 


The scaling of the measurement noises with the associated 
RMS response is a common practice. It is further motivated 
by the form of a human's "observation" noise in man-machine 
studies. The noise -to -signal ratio p is an overall scale 
factor; thus, its actual value is not of large concern. How- 
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ever, a value p = 0.01 tt has been selected from previous 

experience in man-machine systems to represent a nominal 

2 

noise level.* The variances a j.. are picked at the nominal 
point and held constant at these values throughout the 
analysis . 


In any specific design it may or may not be necessary for all 
responses r , , , , , r^ to be insensitive to parameter variations. 
Thus, it is desirable to include within only those responses 
appropriate for desensitive design. The s.^ easily serve this 
purpose through setting = o for these responses. Thus, 


s. 

1 


r (per 



-1 


for desensitive r. 

1 



otherwise 


The weighted trace of the (average) information matrix is 



where W is a diagonal p x p weighting matrix, W = diag (w^). 
The weightings w^ are selected to normalize the information 


*The case where the are not all equal could be considered. 
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matrix at the nominal condition, i.e., at = 1, and to keep 
the analysis in terms of relative sensitivity (with respect to 
the nominal). Thus, it has been found convenient to select 


w. = Y- • 

X '1 


mm ( M ) . . 

os XX 
X 

(M°) 

OS XX 


( 76 ) 


where M is the information matrix evaluated using the nominal 
00 

matrices F , H , D and the nominal feedback gains K . The 
o’ o* o o 

additional scaling factors can be selected to reflect 


1, The relative importance of to the design problem 
as noted through experience or experimentation, or 

2. The relative probability or frequency of occurrence 
of variations in the parameter 


In the present analysis we set 
Yi = 1 i = 1, . . . , p 

to indicate that all parameters are equally important and equally 
likely to vary. 


Thus, unlike some other methods of desensitive controller 
design, the choice of weighting parameters S and W is fairly 
straightforward. In the next section, 3 ^ is appended to the 
original performance cost functional J, given in Equation (70). 
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Design Method 


When the system parameters are subject to variation, a design method 
based upon optimal control theory involves the following two steps, either 
separately or in combination: 

1. Selection of a design point, i.e,, a parameter set at 
which the design is done. Note that need not be the 
nominal point 

2. Selection of a set of optimal gains K at the design point. 

Constraints that are imposed by the physical system, such as limited 
control effort, maximum allowable deviations, etc., are to be satisfied 
in the selection process. 


It would be desirable to achieve both objectives through solving one opti- 
mization problem with a generalized quadratic cost functional. The cost 
functional should reflect the dual goals of performance and de sensitivity, 
and so an intuitive choice is 


J = 




+ BjJ. 


= E{r' Q £} + P2 


(77) 


This is a weighted combination of performance and sensitivity "costs. " 
The selection of both and K could be accomplished by solving a mini- 
max problem, viz, 

J* = min max J(a, K) 

K or “ 
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However, this represents a problem of immense difficulty. An alternate 
approach has been suggested, based on the concept of ’’maximum difficulty”: 


1. Determine the design point on the basis of a maximum 

difficulty criterion. This criterion is dependent on open-loop 
system properties, and so can be found independent of K. 


2. Determine K by solving an ordinary minimization problem, 
min J(ff K). 


This two-step procedure is feasible from a computational viewpoint. The 
first step is discussed in Appendix B. The second problem, finding K, 
is the subject of this effort. But since has not been selected, the 
optimization with respect to K will consider as fixed, but arbitrary. 

As a starting point we will pick “ “ q “ nominal point. Thus, we seek 

K* = argminJ(o;,, K) 

— d 


Another interpretation of the above cost functional is that seeks a K to 
minimize performance at the design point The second term J ^ seeks 

to enlarge the region about in which the gains K remain useful. The 
relative weightings and pg have p^ =1 and Pg chosen so that (after 
finding K*) the resulting does not exceed its minimum value by more 
than a preselected factor 1 + e* where 


1 + e 


K*) 


(73) 
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Thus, we trade off a fraction e > o of performance cost for desensitivity. 
Note that this trade-off need not necessarily be on J^; it could be on 
E{u.^}. E{u.^}, etc. 

Consideration has been given to including a third term 
J_ = E{(r - r )'P(r - r^)} 

within the cost functional J. This would tend to minimize the deviations in 
responses from the original design. Also, it would add terms of the form 
(K - thereby placing constraints on the control gains and indirectly on 

the control effort. It is similar to the uncertainty weighting design of 
Reference 1 , but in a more meaningful closed-loop context. The equations 
that result from appending P 3 J 3 have been developed in detail. But since 
they are more complex than those for J and J„ alone, they will not be 
included in the following sections; in the following sections, we assume 



PROBLEM SOLUTION 


In this subsection a closed-form expression is obtained for the gradient 
matrix ^J/^K that will be used in the subsequent numerical optimization. 

An Expression for the Cost Functional 


The cost functional J, Equation (77), can be rewritten to combine the two 

terms J and J„, The result is 
1 ^ 
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I 


= Efx ' Q X 

XX — 


P 

+ S 
i=l 




Q 


O’; O; 


H.i 


„ P 
+ 2 2 

i =1 


Q a .] 

xa^ - i-’ 


(79) 


•where 

5x 

C7 . = — ^ 

-1 bQ?£ 


is the closed-loop sensitivity vector at the point a = Defining 

F = F + K (80) 

as the actual closed-loop system matrix, o;-^ is seen to satisfy the 
differential equation 



i = 1, 


.,P 


(81) 


where we recall that G 


1 


and G 2 are not functions of a. 


The weighting matrices in Equation (79) are 


Q 

XX 

p 

= pTH+DK)'Q(H+DK) + 2 (H.+D.K) ' S.(H.+D.K) 

(82a) 

1=1 


Q 

(TiOi 

= p_(H-l-DK)'S.(H+DK) 

(82b) 

Q 

xa; 

= P_(H.+D.K)' S.(H+DK) 

^ J. 1 1. 

(82c) 


where the matrices 

S. = w. S i = l, ...,p (83) 

11 » It' 


The state vector, x, and sensitivity vectors can be combined into an 

(n + 1) • n augmented system. Defining 
P X 
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(84) 


Xa -1’ 

as the augmented state, the resulting system equation is 
2iA ^2A 

where 






F 0 . . . . 

0 


^2 

F^ F 



0 

F 0 F 


; G = 

• 

• ^ 
• 


A 


F 0 .... 

0 F 


0 

^ P 





It is important to point out that this augmentation is done purely for 
analytic simplicity. Fortunately, it will not be necessary to solve large 
dimensional linear (or nonlinear) matrix equations in the ensuing opti- 
mization process. 

In terms of X^, the cost functional J can be written as 
where 
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Note that there are no cross -terms between cj- and a.. We define 

1 .1 



as the covariance matrix of the augmented state. 


<=A =■ 

Cov[x^] = E{x^Xa^ 

(86) 

a full matrix with 


C 

XX 

E{x x'} 

(87a) 

C 

E{X CTj^} 

(87b) 

C 


(87c) 

CTiO. 
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Using the cyclic property of the trace in Equation (85) and substituting 
Equation (86) yields 

J = tr [Q^C^] (88) 


From Equation (84) it is seen that satisfies the Lyapunov equation 




(89) 


or, equivalently. 


~ J*o ® ^ *^2A^2A ® (90) 

Incorporating Equation (90) into Equation (88) and manipulating terms gives 

J = °2A °2a1 <«' 


where 

L.=| e AQ.e A dt 
A o A 


satisfies the linear equation 

^'a ^a^^a^a +'3a ’ 0 


(92) 


Equation (91) is better to use in the analysis than is Equation (88) since both 
and in Equation (88) depend on K, whereas only L^ in Equation (91) 
depends on K. 


Equations (91) and (92) may be further simplified by investigating the special 
block structure and form of L. : 
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L L L 

XX xa. Xo-r 

1 




Xar 


0 

L 


cTocr 


2^2 


L 0 


. L 


XCT 


. 0 


cr CT 
P P 


It is easy to show, by expanding Equation (92), that the cross-terms 
L =0. The equations for the nonzero components are 

UiO. 


0=FL +L F+Q 


CT-Oi O-Oi 


a.Oi 


i = 1, . . . ,p 


0=FL +L F+Q +F.L i = l,,..,p 

Xct^ xcr^ 1 


(93a) 

(93b) 


_ P 

0 = FL +L F + Q + s (L F. +F/L' ) 

XX XX XX . ^ xa. 1 i xq. 

1=1 1 1 


(93c) 


These equations are solved in the order given. We first solve for L ,* 

CTiCTi 

then for , and then for L . The matrix L is all that is 

XX XX 

needed in evaluating J because of the sparse form of ^2A*^2A' from 

Equation (91), 


wi 

+Note L = — r L so that only the equation for a- need be solved. 

CTiCTi W] CTjCTj i 
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J = tr a^G'1 


(94) 


The sequence of lower order matrix equations to compute the covariance 
terms (62) and (63) is given by 

0 = FC +C F" +G„G' (95a) 

XX XX 2 2 

0 = FC +C F'+C F/ i = l,...,p (95b) 

X(T. X(J. XX 1 

1 

0 = F C +C F"+F.C +C" F/ (95c) 

CTiCTj a^O-j ^ 3 

for i = 1, . , .,p; j = i, . . .,p 

Thus, on the basis of Equations (93) and (94), only (p + 2) linear equations 

of order n^ need be solved to obtain the value of J at a given K. Moreover, 

all of these equations involve the same matrix F; the covariance equations 

(95) all involve F'. It is shown below that a modified Bartels -Stewart 
5 

algorithm can be used to reduce greatly the computational burden asso- 
ciated with evaluating Equations (93) through (95). 

The optimization problem is therefore to find the constant gain matrix K* 

to minimize tr TL G„GiT (at the nominal point a , = a ), such that 
^ XX 2 2 ^ — d — o 

Jl(£o» K*) ^ (1 + e) Ji(a^, K^) 

Note that this approach constrains the control u(t) to be of the form 
u(t) = K X (t) 

There is no feedback of the sensitivity vectors as in the sensitivity 
vector approach. 
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Gradient Expressions 


It is not possible to find a closed-form expression for the gains that mini- 
mize J, except when ^2 “ For this reason, the optimization will be 
approached via some form of gradient algorithm. Closed-form expressions 
for the gradient matrix fjJ/gjK will be of great advantage in this process, 
since numerical evaluation of gradients would be extremely time-consuming. 


The gradient 7 t^J = bJ/jjK is evaluated using a technique of Kleinman for 
derivatives of trace functions. From (93c), the first order variation in 
to a change K -* K + §K is (note that F is a function of K), 


0 = F '6L. 


XX 


6L F 

XX 


6K'G/L 

1 XX 


L G, 6K 

XX 1 


+ 


XX 


P 

S 

i=l 


(6L F. 
xa. 1 
1 


F/ 

1 


6L' ) 


(96a) 


where 

0 = F' 6L + 6L F + 6K'G; L + L G, fiK + sQ 

xct. xct. 1 xct. Xa. 1 Xct. 

11 11 ^1 

+ F/ 6L 
1 CT. a,. 


(96b) 


and. 


0 = F'6L 

for i = 1, . . . , p. 
given by 


+ 6L F + 6K'G/L +L G, 6K + 6Q (96c) 

a.a. 1 a.a. a-a. 1 a-a- 

The first order variations in the components of are 


79 



*^xx ^ p^[«K'D'Q(H + DK) + (H + DK)'QD8K] 

+ p2 ^ 'S.D. 6K] 

i=l 

6Q = Bo r(H. + D.K)'S.D6K + 6K^D/S.(H + DK)] ( 

xa. /‘"ill 11 

1 

&Q = So [SK'D'S.(H + DK) + (H + DK) 'S. DfiK] ( 

CT^CTj^ I 1 

The matrix F is the closed-loop system matrix which is required to be 
stable for any choice of K. Equations (96a) through (96c) can therefore 
be written as equivalent integral expressions, 

aL = r“ e^'Yj-jK'G/L +L GaK+ 6Q ]e^^dv 

a.<j. -o 1 o.a. a.a. 1 a. a- 


aL = r”e faK'G/L +L G, aK + aQ +F/aL ] e dt 
Xct. *^0 1 xa. Xct. 1 Xcr. i a-a. 


1 1 


where a similar expression for SL^^ can be written directly from 
Equation (96a). 


Since aJ = tr.[eL GoGi], we substitute for aL its integral expression, 

XX Z ^ XX 

Substituting further the integral expressions for aL and aL , using 

‘^i'^i 6 

the cyclic properties of the trace and Kleinman's lemma, one obtains, 
after tedious manipulation. 


V^J = 2G/L C + 2b,D'Q(H + DK)C 
K 1 XX XX 3 


+ 2G/ rCL C' +L'C +L C ] 
1 xa^ xa. xa. a^aj 
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p 

+ 2 Bo S rD/S.(H + DK)C ‘ + D 'S.(H. + D.K)C 

^2. , ^ 1 1 X(j. 11 1 Xct. 

1=1 1 1 

+ D'S.(H+DK)C + D/S.(H. + D.K)C I (98) 

1 a- a. Ill 1 xx-^ 

1 1 

Equation (98) will serve as the basis of a gradient algorithm to minimize J, 
The computational requirements to compute and J at a given K are 

now the major issues. 

COMPUTATIONAL ALGORITHMS 

This section describes the numerical schemes for computing J, 

M for a given feedback. 

00 

General Overview 


In order to numerically evaluate the cost functional J we must solve (p + 2) 
linear matrix equations (93a) through (93c). To evaluate the gradient Tr^I, 
an additional (2p + 1) equations need be solved for C , C , and C 

XX X(7j^ 

as in Equations (95a) through (95c). Note that the cross - 
correlations C i ^ j are not needed to compute . Examination of 

CTj^CTJ 

these 3p + 3 equations reveals the following: 

1. All equations involve F or F.' Half of the linear equations 
are adjoint to the other half. In particular. Equation (95) 
is adjoint to Equation (93). 
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2. A total of p + 3 equations have symmetric solutions, whereas 
the other 2p equations for L and C have nonsymmetric 

XO'j^ 

solutions. 

A significant reduction in computation time can be brought about by several 
modifications to the linear equation algorithm of Bartels and Stewart, This 
algorithm is well suited to the efficient sequential solution of 

A'X + XA = C. (99) 

with different right-hand sides The existing algorithms and available 

programs are geared for symmetric and hence, symmetric X. By a 

slight modification a skew -symmetric can be handled as well. Thus, 

the two major objectives for sensitivity design are 

1, To solve Equation (99) when (and therefore, X) is non- 
symmetric, and still take advantage of the saving in computer 
time and efficiency afforded by a symmetric problem, and 

2. To solve the adjoint equation 

AX + XA' = a (100) 

using a computer program written to solve Equation (99), 
Modifications to the Bartels -Stewart Algorithm 

The Bartels -Stewart algorithm for solving Equation (99) in the symmetric 
(or skew-symmetric) case first reduces A to an upper Schur from via an 
orthogonal transformation Q. The resulting matidx A is of the form 
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rw 

'^ll 

^12 • 

. . A, 

Ip 

A = Q' A Q = 

0 

A 

•^22 • 

■ ■ *2p 


0 


/XV 

' . A 

PP 





where each submatrix A.. 

u 

is at most a 2 x 2, 


/X/ 

<=11 

^12 • 


C. = Q'C. Q = 
1 1 

r>* 

<=21 

• 

°22 • 

• ■ "2P 

• 


m 

_°pi 


• 

• 

c 

pp 

and 





*11 

*12 • 

. . X, 

Ip 

X = X Q = 

*21 

/-w* 

V 

22 • 

• • *2P 


X , ... X 

pi pp 

then Equation (99) is equivalent to 

n-» 

A X + X A = C. 

1 


( 101 ) 


Since the partitions of C. and X are conformal with A, expanding Equation 
(101) gives 




k-1 


X-1 


Aj^ X,. . + X,_ A . . = C 


- 2 A,'. X. - 2 X, A. 

kX kx ZJi ki ^3 k. li 

X ”* 1 1 j Pi k "" • i P 


( 102 ) 
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These equations can be solved sequentially for X , X_ , . ,,X „ 

11 ^1 ^Ct pi 

. . . ,X ; we then fill in the upper part of X by symmetry 
PP 






or skew -symmetry 





as the case may be. The ''mini-systems” of Equation (102) are solved via 

r>u 

a simple algebraic equation program. Thus, with X calculated, the 
solution X to Equation (99) is 

X = Q X Q' 


Once the real Schur form of A and Q have been calculated, they may be 
saved and reused to solve the same equation (99) with different Cj^, They 
are also used in the iterative refinement of the computed solution. This 
is the forte of the Bartels -Stewart algorithm. 

Solution of the Nonsymmetric Case 


Any general C. can be written as a sum of a symmetric part C. 

^ ( 2 ) ^ 

skew symmetric part where 


( 1 ) 


and a 
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= i(C. - CO 

For a general matrix the linear equation (99) can be solved by adding 
the solutions of the two equations 

A = (103a) 

A'X^^^ + X^^^ A = (103b) 

Thus, X = X^^^ + X^^^, where we note that X^^^ is symmetric and X^^^, is 
skew -symmetric. 

Once an upper Schur form of A is available, say from a previous solution of 
Equation (99), we can call the computer program twice and solve for a non- 
symmetric solution. The time required to solve the linear equation (99) 

i-W 

once the matrices A and Q have been found is about 40 percent of that 
required to solve the equation for the first time. Therefore, a nonsym- 
metric solution is obtained in about 80 percent of the time needed to solve 
the equation once. 

Solution of the Adjoint Case 

The adjoint equation (100) could be solved rapidly if somehow we could 
obtain A' in upper Schur form from A. Consider Equation (100) where 
Q and the upper Schur form are given. Pre- and post -multiply this 
equation by Q ' and Q, and note that QQ ^ = 1 gives 

Q ' A Q (Q 'X Q) + (Q "X Q) Q 'A'Q = Q 'C. Q 
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or 


A X + X A' = C. 


(104) 


where A' is in lower Schur form. Now we need A' in upper Schur form 
to use a program written to solve Equation (99). The matrix A* can be 
transformed to upper Schur form by applying a symmetric, orthogonal 
transformation T to Equation (104) where 


T = 


0 


0 


Pre- and post-multiplying Equation (104) by T gives 
TAT(TXT) + (TXT)TA'T = TC.T 


or 


A/X + X A = TC.T 
11 11 1 


(105) 


r-^ 

where A^^ is now in upper Schur form as required. A summary of the 
steps needed to solve Equation (100) is as follows: 

Transpose A . 

rsj 

Obtain A^ = T A ^ T, so A^^ is A ^ with its rows and columns 
written in reverse order. 


1 . 

2 , 

3, 

4, 

5, 


Obtain T a T. 


Solve for X^^ using the same algorithm as for Equation (99), 


Obtain X =T(QX^Q0T. 


86 



A general-purpose computer program, AXPTA, has been written to solve 
the linear matrix equation via the Bartels -Stewart approach. It has the 
features to solve the symmetric, skew-symmetric, and adjoint cases 

r*./ 

taking advantage of previously obtained A and Q. 

Computational Requirements for J, and 


As noted above, p + 3 symmetric and 2p nonsymmetric solutions are needed 
in evaluating J and Taking advantage of the algorithm modifications, 

these symmetric equations can be solved in an "equivalent" computational 
time of l+0.4(p + 2) = 0.4p + 1.8 linear equations. The remaining 2p 
equations with nonsymmetric right-hand sides can be solved in the equiv- 
alent of 2p X 0,8 = 1, 6p linear equations. The total computation time is 
then « 2. Op + 1. 8 linear equations. Thus, to obtain J andv^ J for p = 3, 

XV 

we need solve the equivalent of roughly eight linear equations. This is 
comparable to the time required to solve one (n -dimensional) Riccati 
equation. Thus, the computations of J and its gradient (at each iteration) 
are not excessive. 

The unit time information matrix M is not required explicitly in the opti- 

CO 

mization algorithm. However, it is useful to monitor M and the disper- 

00 

Sion matrix 



as the algorithm proceeds through the iterations to see how the uncertainty 
regions for increase. Certainly, one would wish to compare or 
at the optimal point K* with their initial values at K^. In addition, the 
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(initial) diagonal elements of M are used in forming the weighting factors 

00 

Wj^ in Equation (52). 


The p X p information matrix M is 




b£ 5£ 

boi^ 


The response sensitivities are 

b£ 

— = (H. + D.K) X + (H + DK) a . 

1 1 — —1 

Substituting these into Equation (107) gives for the ij element of M , 

00 

(M ).. = tr r(H + DK) ' S (H + DK) C + (H. + D.K)'S(H. + D.K)C 

00 13 a.CT. 11 3 3 XX 

o 

+ (H + DK) ' S (H. + D.K)C + (H. + D.K) " S(H + DK)C ' ] 

3 3 xa^ 11 xa^-* 


(107) 


Notice that C , C and C will have already been evaluated while 
XX* Xcj. cr.CT. 

1 11 

computing Therefore, it is only necessary to determine C for 

K CT£aj 

i ^ 3 from the Equation (95c). Only the terms for 3 >i need be 
computed because of symmetry. This requires solving an additional 
p(p - l)/2 = 3 linear equations, all with the same system matrix F. 


Iterative Algorithm for Finding K* 

Returning to the optimization problem for K*, two gradient algorithms are 
proposed in this section. The first algorithm was tried first, primarily 
to check the validity of the overall approach to the sensitivity problem. The 
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second method is based on the conjugate gradient algorithm of Fletcher- 
Reeves. As noted earlier, the computer programs solve for K* at a 
given design point assumed here to be the same as 

Successive Substitutions Scheme 


This is a simple iterative scheme that has no proven convergence prop- 
erties. It is motivated by old algorithms for the optimal output feedback 
problem. The idea is to find a set of gains ^ that would result in 
= 0. The choice is based on matrices computed at iteration x. The 
algorithm is as follows: 

1. Set < 4=0 select an arbitrary initial (stabilizing) gain 

Usually = K^. Set = ». 

( ;) ^ ( i;) 

2. Compute L , L , L and J(K' ^') = j' 

0^.0** Xo" • XX 

11 1 

3. Check stopping condition ^^1 < TOL 

If satisfied, stop. 

( a) 

4. Compute C , C , C and the gradient V . J 

XX’ XCT.’ CT.a. ^ j 

( L) 

5. Compute the gain increment AK that would make 

^ = 0 assuming all other matrices remained con- 

stant. As seen from Equation (98) this is a near-impossible 

task. Thus, the gains K in the summation term are set to 
(») 

K ; i.e., they are fixed, and we find 
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6. 


Select new gains ^ + b where b = ^ 

initially. If the cost ^ does not decrease, a smaller 
step is taken by reducing b. 

7. Set X = jt + 1 a-nd return to Step 2. 

The above algorithm is essentially a successive substitution scheme for 
solving an equation of the form x = f(x). Such a scheme is convergent 
only if the slope of f is < 1. Thus, the given algorithm is expected to 
converge when the optimum gain K* is close to the initial gain 
Unfortunately, the convergence rate of this method was found to be very 
slow, with considerable oscillation in K when near the optimum. Con- 
vergence, to less than 1 percent error, was usually attained in 10 to 
18 iterations. 

Conjugate Gradient Method 


A conjugate gradient scheme was picked as an alternative to the above 
method. The steps are as outlined: 

1. Set 1 = 0 select initial gain Usually, we set 

= K^. Pick M = p as recycle index. 

2, Compute 
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3. 


If J 


(j») («) 

Compute the gradient V J at the current gain K * 

is too small, then stop. 


4, Compute the current descent direction 

“ 

where a° = 0, and on subsequent iterations 
a^ = tr {[ 

tr {[ ~ ' s^-^ " 

Normalize s^"^^ so that || || =1. 


5. 


Compute the current step size 


using a one -dimensional search 


b"^ = arg 


mm 

b 


J + b 




This is done via a quadratic interpolation scheme. 


6. Compute new gains 


+ 1) = j^(x) + 


and the cost J 


U + 1 ) 


7. Check convergence tests. If passed, then stop. 

8. If jt < M set X = jt + 1, Otherwise, set l = o 

Go to Step 3. 
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The convergence rate of the conjugate gradient method has proved 
superior to the first noted algorithm.* Typically, three to eight iterations 
have been needed for convergence, CPU time on 360/65 of about 6 to 7 min 
utes. The critical parameter for convergence is the initial guess for b in 
the one -dimensional search (Step 6). The present ad hoc guess is 

, II II Min Abs element in 

|3 ^ • / -\ 

Max Abs element in K' ^ 

PRELIMINARY EVALUATIONS 

The 15 -state residualized model of the C-5A longitudinal dynamics as 
described in Reference 1 was used for evaluation purposes. The response 
vector used for the design of the nominal controller is ninth order and is 
defined by 

r^ = Bending moment at wing root 

r„ = Torsion moment at wing root 

r^ = Bending moment rate 

r • = Torsion rate 
4 

r^ = Aileron displacement 

r_ = Inboard elevator displacement 
b 


It) 

*A combination of the two algorithms to give was tried. However, 

it gave mixed results and so was not pursued. 
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= Modified aileron rate 
r = Inboard elevator rate 

O 

rg = Control follower response 

The weighting matrix Q used in the quadratic synthesis design is diagonal 
with 

Q ={0. IE - 9, 0. IE - 9, 0. IE - 10, 0.32E + 8, 0.0, 0. IE + 7, 

0. IE + 5, 0. 2E + 6} 

The parameter p vector a for the C-5A example is 
Of ~ u)£* ^^1/7 ^ 

where 

= dynamic pressure uncertainty factor 

(j)^ = structural frequency uncertainty factor 

= stability derivative uncertainty factor 
f 

The design process is conducted at the nominal condition, a.^ = 1, The 
gradient matrixes for i = 1, 2, 3 were computed numerically. 

They are assumed to be constant within the range of parameter variation 

The sensitivity reduction problem formulated via the Information Matrix 
approach leaves very few free parameters to be chosen. This is by design, 
since we have elected to minimize the free parameters in order to mini- 
mize the number of design iterations. As a result, one needs only select 
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1. Which of the responses are to be desensitized, i.e., 
which responses n(t) have Sj^ ^ 0, and 

2. The value of Pg* giving a relative weighting to the sensitivity 
cost versus the performance cost. 


As a first step in evaluating the technique, the effect of varying on the 
bending moment response, r^, was investigated. In other words, 
s ^ = 1, = 0, j = 2, 9. Wj^'s were set to unity to reflect equal impor- 

tance on all uncertain parameters. Note that the weighted trace of the 
information matrix in this case is 



As a means of evaluating the effect of varying ^ measure of local 
sensitivity was defined as 



where the term r^ ^ refers to the bending moment response at the nominal 
condition. This is slightly in variance to the more classical definition of 
sensitivity which is given by 


bE(r^) 


“i 


E{r^} ^“i 
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1 

E£r^5 



2pE{r^)t 

E{r^) 




where n is the correlation coefficient between r and at the nominal 

bar 

point. Hence, the measure used is proportional to the more classical 
definition. 


Figure 25 presents a plot of S ^ versus As can readily be seen, there 

Of 4 

is a reduction in sensitivity, particularly with respect to structural fre- 
quency uncertainties . 


The effects of varying Pg were also ascertained with respect to J^, the 
performance cost index, control activity measured in terms of aileron dis- 
placement, and the identifiability of the uncertain parameters measured in 
terms of standard deviations of the estimates. These results are plotted 
in Figures 26, 27, and 28, respectively. As expected, the results are as 
follows: 

• Performance as measured by the cost index increases 
with 62 * 

• Decreased sensitivity requires increased control activity. 


Uncertain parameters are harder to identify as Pg 
increased. 
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© DYNAMIC PRESSURE q^(o^) 

© STRUCTURED FREQUENCY 
© STABILITY DERIVATIVE M^^ (og) 


■-0 ^ 


— o 

-CL © 

*" 0 — _r> 





-O 







/ 




o DYNAMIC PRESSURE q^(a^) 

© STRUCTURAL FREQUENCY W^(a2) 
STABILITY DERIVATIVE M^^(a2) 


1 


1 


1 
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Since the sensitivity measure employed is only valid over small variations, 
a linearity test was also performed. The parameters were varied 
independently from -0. 2 to 0. 2 in steps of 0. 05. The normalized incre- 


AcT: 


mental rms response 






Aof. 


and linear "prediction" 




1,0 


1. O 


are listed in Table 7, The results are shown for parameters 1 and 3 only. 


For oTo, the analysis showed that the correlation between r. and 1 / 

2 1 Xba 


was very small. Thus, although E' 



is large, the effect of 


variations in «2 actual response is small. As seen in Table 7, the 

changes are approximately linear. 


TABLE 7. INCREMENTAL VERSUS LINEARIZED PARTIALS COMPARISON 


1 

Parameter 

Variation 

Actual 

Performance Change 

"Predicted" Change 

Aofj^ 

Act 

^1 

""i,© 


bar. 

AOTi y 
/ ^1,0 


i=l 

i=3 

i=l 

i=3 

0. 05 

0.031 

0.015 

0.031 

0, 0142 

0. 10 

0.061 

0, 031 

0,059 

0,028 

0. 15 

0.090 

0.049 

0,086 

0,042 

0.20 

0. 1193 

0. 068 

0,1111 

0,056 

-0.05 

-0,0327 

-0.0137 

-0,0322 

-0.0144 

-0. 10 

-0. 0679 

-0. 0262 

-0, 0654 

-0. 0291 

-0, 15 

-0.1055 

-0.038 

-0, 1009 

-0. 044 

-0.20 

-0.1460 

-0. 0483 

-0. 1357 

-0,0591 
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With the effects of varying established, a refined design was under- 
taken and evaluated. This is discussed in the following section. 

C-5A CONTROLLER DESIGN EVALUATION 

As discussed in the preceding section, the Information Matrix design 
technique was formulated to limit the number of free parameters which 
the control system designer must manipulate. As presented here, the 
designer has freedom to vary the scalar term Pg which weights the sen- 
sitivity reduction. The designer also may select which of the system 
responses he wishes to desensitize by manipulation of the binary variable 
s.. The effect of varying p„ on one system response was discussed in 
the previous section. This section will discuss the effect of varying both 
^2 and Sj^ with the purpose of obtaining an insensitive C-5A control 
system. 

Design Approach 


In order to limit the freedom on selection of additional constraints 

were imposed; 

1. The value of of the insensitive controller must be less 
than 1,2 times the of the nominal controller. 

2. Aileron and elevator control activity for the insensitive 
controller must be less than two times the controller activity 
of the nominal. 
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With these constraints imposed, computer runs were made to study the 
effects of 

1, Sensitivity weights on rates B^, to aid in desensitizing 
the bending and torsion moment responses, 

2, Weighting aileron control responses r^ to keep control 
effort from rising too rapidly, 

3. Weighting the controller follower response r^ as an attempt 

to desensitize uj , r , and 
sp’ ^sp’ 

4. Various choices of Pg- 

A total of 15 different cases were studied, including the nominal. The 
cases, described in Table 8, may be grouped into five categories: 

1. Variations in Case 2A, 2B, 2C; Case 3, 3A, 4; Case 7,8 

• « 

2. Sensitivity weights on B^, T^: Case 7,8 

3. Sensitivity weight on 6a: Case 2A-2F, 3, 3A, 4 

4. Sensitivity weights on 6^^: Case 2A-2F, 6 

5. Variations in the (pseudo) noise/signal ratio cp on 6a to 

study more closely effect of control: Case 2C-2F. 
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TABLE 8, PARAMETER ESTIMATE STANDARD DEVIATIONS AND OTHER 
PERTINENT INFORMATION FOR DIFFERENT GAINS 


Run If 

02 

Diagonal 
elements 
of M„ 
with 
nominal 
gain K° 

Weights 

on 

trace 
of Mco 

Initial 

estimates 

Lower Bounds On 

Ratios of 
responses 
with K* and K° 

Any 
weights 
on Tg 
and 
i.e 
and Y 

Any 

weights 

on 

control 
response 
i.e., rg 

Any 

weight 

on 

control 

followed 

response 

Signal/ , 

noise ' 

ratio on 
control, r^ 

(% pseudo 
noise) 

Final 

estimates 

Volume of 

Volume of 
uncertainty 
ellipsoid at 
K*x lO"^ 

Cost, J 

1 

( Ratio) 

of S-D of 
parameters 
with K° 

of S-D of 
parameters 
with K’^ 

uncertainty 
ellipsoid at 
K° X lO"^ 


^2 



1 

NA 

NA 

NA 

NA 

NA 

NA 

NA 

88. 7 

1.0 

1.0 

NA 

NA 

NA 

NA 

(Nominal) 








(1.0) 









116.2 

0.40 

0.0933 

0. 1 146 










2A 

0.4 

538.2 

0.08 

0. 0431 

0.0499 

3.36 

8.035 

91.0 

0.977 

0.944 

No 

Yes 

Yes 

31.8 



49. 1 

0.95 

0. 1436 

0. 1551 



(1.026) 






(1.0) 



120.8 

0.42 

0.0913 

0. 1269 










2B 

0.6 

559.6 

0.09 

0.0423 

0.0590 

4.82 

12.518 

93.6 

0.940 

0.939 

No 

Yes 

Yes 

31.8 



50.8 

1.00 

0.1409 

0. 1512 



(1.055) 






(1.0) 



120.8 

0.42 

0.0913 

0. 1110 










2C 

0.5 

559.6 

0,09 

0.0423 

0.0473 

4.82 

6.770 

93.6 

0.952 

0.924 

No 

Yes 

Yes 

31.8 



50. 8 

1.00 

0. 1409 

0. 1571 



(1.055) 






(1.0) 



78.0 

0.64 

0.1155 ' 

0. 1229 










2D 

0.5 

253.6 

0.20 

0.0634 

0.0670 

10.52 ' 

16.434 

92. 1 

0.940 

0.943 

No 

Yes 

Yes 

3. 18 



40.8 

1.00 

0. 1432 

0. 1569 

1 


(1.038) 






(10.0) 



87.0 

0.558 

0.1086 

0. 1177 










2E 

0.5 

325.2 

0. 149 

0.0557 

0.0583 

7.47 

11.21 

90.38 

0.954 

: 0.955 

No 

Yes 

Yes 

10.6 



48.6 

1.0 

0. 1448 

0. 1548 



(1.020) 






(3.0) 



94.9 

0.558 

0. 1036 

0. 1298 










2F 

0.5 

381.8 

0. 149 

0.0513 

0.0630 

5.79 

15.937 ' 

91.72 

0.947 

0.949 

No 

Yes 

Yes 

13.2 



48.8 

1.0 

0. 1443 

0. 1593 



(1.034) 






(2.0) 



104.9 

0.40 

0.0989 

0. 1372 










3 

0.6 

497. 7 

0.09 

0.0449 

0.0568 

4.90 

23.49 

95.4 

0.937 

0.910 

No 

Yes 

No 

31.8 



41.2 

0.95 

0. 1579 

0. 1968 



(1.076) 






(l.O) 



104.9 

0.40 

0.0989 

0. 1448 










3A 

0. 7 

497. 7 

0.09 

0.0449 

0.0688 

4.90 

38.90 

94.47 

0.9126 

0.9026 

No 

_ Yes 

No 

31.8 



41.2 

0.95 1 

0, 1579 

0. 1995 

1 



(1.065) 






(1.0) 
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TABLE 8. PARAMETER ESTIMATE STANDARD DEVIATIONS AND OTHER 
PERTINENT INFORMATION FOR DIFFERENT GAINS (concluded) 


Run II 

02 

Diagonal 
elements 
of M„ 
with 
nominal 
gain K° 

Weights 
on 
trace 
of M« 

Initial 
estimates 
of S-D of 
parameters 
with K° 

Lower Bounds On 

Ratios of 
responses 
with K- and K° 

Any 
weights 
on rg 
and 

i.e., Bj^ 
and Y 

Any 

weights 

on 

control 
response 
i.e., r^ 

Any 

weight 

on 

control 

followed 

response 

Signal/ 
noise 
ratio on 
control, T 5 
(% pseudo 
noise) 

Final 
estimates 
of S-D of 
parameters 
with K- 

Volume of 
uncertainty 
ellipsoid at 
K°x lo"^ 

Volume of 
uncertainty 
ellipsoid at 

K- X 10? 

Cost, 

( Ratio) 


T2 


r^ K° 



104. 1 

0.40 

0. 0989 

0. 1296 










4 

0.4 

497.7 

0.08 

0.0499 

0.0546 

4.90 

17.43 

91.4 

0.961 

0.931 

No 

Yes 

No 

31.8 



41.2 

0.95 

0. 1597 

0. 1867 



(1.03) 






(1.0) 



59.8 

0.68 

0. 1385 

0. 1452 










5 

0.2 

169.3 

0.20 

0. 0796 

0.0921 

15.03 

50. 78 

90. 1 

0.931 

0.967 

No 

No 

No 

NA 



40.3 

0. 80 

0. 1632 

0.1780 



(1.019) 









67. 1 

0.68 

0. 1258 

0. 1277 










6 

0.2 

290.6 

0.20 

0.0603 

0.0824 

13. 14 

25.54 

90.0 

0.960 

0.958 

No 

No 

Yes 

NA 



41.5 

0.80 

0. 1555 

0. 1573 



(1.018) 









127.2 ' 

0 800 

0.0942 

0. 1008 










7 

0.2 

1432.8 

0.043 

0.0281 

0.0331 

0 689 

1 . 259 

91.2 

0.920 

0.955 

Yes 

No ' 

No 

NA 



101.6 

0. 700 

0.0992 

0. 1065 



(1.027) 









127.2 

0.800 

0.0942 

0.0978 










8 

0. 1 

1432.8 

0.043 

0.0281 

0.0306 

0.689 

0. 960 

89.5 

0.945 

0.978 

Yes 

No 

No 

NA 



101.6 

0. 700 

0.0992 

0. 1036 



(1.0009) 









59.8 

0.612 

0. 1385 

0.1484 1 










9 

0.3 

169.3 

0.237 

0. 0796 

0. 0994 

15.03 

70.0 

91 3 

0.922 

0.941 

No 

No 

No 

NA 



40.3 

0.800 

0. 1632 

0. 1864 



(1.029) 










The performance of the resulting 14 controllers plus the nominal controller 
are presented in Table 9. The 14 controllers were then compared to deter- 
mine which one would be evaluated with the design criteria defined in 
Reference 1. The Case 3A controller was chosen based on conditions that 

1. All the design specs were satisfied plus the additional 
imposed constraints on and control activity, and 

2. The identifiability of the three uncertain parameters 

(q^, Wj, was reduced the most. This was measured 

by the volume of the uncertainty ellipse (which is approxi- 
mately the determinant of the dispersion matrix) and the 
standard deviations of the uncertain parameters. 

The gains for Case 3 A are given in Table 10. With respect to the other 
variations that were investigated, it was found that 

1. Weighting the control follower response offers no advantage. 

This is to be expected since the control follower response 
is only valid at the nominal condition. Its purpose is to 
achieve a specific control configuration (i.e., (,e^ = 0,5q) 
which at the nominal produces desirable short period 
frequency and damping characteristics. 

At other than the nominal condition, the control follower 
response will not produce the described short period 
frequency and damping characteristics. 
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TABLE 10, GAIN MATRIX FOR RUN #3 A 


11 

-1.5756E-05 

12 

5.5 747E-06 

13 

3.2549E-05 

14 

2.0883E-05 

15 

-2.0162E-04 

16 

-7. 04 19 E -04 

17 

-2. 1381E-01 

18 

1. 1257E-03 

19 

2.5313E-05 

110 

-2.4873E-04 

111 

8.9334E-05 

112 

-8.5690E-06 

113 

-1.9357E-04 

114 

1.3001E-04 

115 

-4.5946E-07 

21 

1.9712E-04 

22 

5.0306E-04 

23 

-2.8724E-04 

24 

2.6991E-04 

25 

4.7223E-04 

26 

-6.2514E-03 

27 

3.4040E-03 

28 

-1.5226E-01 

29 

6.7437E-06 

210 

6. 7847E-04 

211 

-1.2334E-03 

212 

-3. 1175E-05 

213 

4.4720E-04 

214 

8.4716E-04 

215 

1.3936E-03 



2. There is no significant change in system response as the 
signal/noise ratio on 6 is varied. 

Si 

• • 

3. The effect of including and T^^ in the sensitivity computa- 
tions appears to have little effect on the resultant and T^^ 
responses , 

Design Evaluation 


The Case 3A Information Matrix controller was evaluated on the 15 -state 
Case 4R residualized C-5A model at the six evaluation conditions. These 
conditions, chosen in the Reference 1 study, are as follows: 


1. 

Nominal condition: 

a = (1.0, 1.0, 1.0) 

2. 

Worst Case 1: 

a = (1.25, 0.75, 0.8) 

3. 

Worst Case 2: 

a = (0.5, 1.0, 1.2) 

4. 

Independent Variation 1: 

a = (1.0, 1.0, 0.8) 

5, 

Independent Variation 2: 

a = (1.0, 0. 75, 1.0) 

6. 

Independent Variation 3: 

Of = (1. 25, 1,0, 1.0) 


The performance of the Information Matrix controller is tabulated in 
Table 11 at the six evaluation conditions. Figures 29 through 33 graph- 
ically portray the tabulated data for each of the design specifications versus 
the performance of the nominal controller. 

The three criteria defined in Reference 1 were used for evaluating the Infor- 
mation Matrix controller. The criteria may be briefly described as follows: 
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TABLE 11 


INFORMATION MATRIX CONTROLLER PERFORMANCE 
EVALUATION MODEL— CASE 4R 


1 Specification 

Description 

Criteria 

Nom 

=1.0 
uif =1.0 

M =1.0 

WCl 

if =1.25 

=0.75 

M =0.8 
^f 

WC2 

=0.5 

UJf =1.0 

M =1. 2 
^f 

PI 

qf =1.0 

U)j^ =1.0 

M =0.8 

Wf 

P2 
=1.0 
=0.75 
M =1.0 

P3 

=1.25 

Wf =1.0 

M =1.0 
^f 

Maneuver Load 
% Change 

B 

< -30% 

-45 . 4% 

-32. 1% 

-73.9% 

-45.7% 

-45.3% 

-32.5% 

Gust Load 


B 

< -30% 

-40.6% 

-29. 1% 

-63.3% 

-37.7% 

-40. 1% 

-32.0% 

Alleviation 










% Change 


T 

<+ 5% 

-37.9% 

-20. 7% 

-59.6% 

-32.9% 

-36 . 3% 

-29.8% 




> 1.5 

2.36 

3.01 

1. 17 

2.26 

2.30 

3,08 ’ 

Handling 


sp 

rad/sec 







Qualities 


r 

0 

1 

O 

oo 

0.72 

0.86 

1 0.677 

0.755 

0. 768 

0.728 



^sp 

sec ^ 







Stability 










Margins 










Gain: 

aileron 

> 6db 

QD 

00 

00 

00 

00 

00 


elevator 


32db 

18db 

=1^ 

32db 

2ldb 

27db 

Phase: 

aileron 

> 45° 

CO 

00 

03 

00 

00 

00 


elevator 


118° 

118° 

125° 

115° 

125° 

120° 

Surface 


6a 







0.00022 

Activity 










RMS 


6a 



0.0012 




0.0010 




NA 







(rad 


6e 



0.0023 




0.0019 

rad/sec) 












6e 



0.0061 




0.0053 










NOMINAL 

INFORMATION 

MATRIX 


Figure 29. Case 4R Maneuver Load Performance 



EVALUATION CONDITIONS 


NOMINAL 

INFORMATION 

MATRIX 


Figure 30 


Case 4R Gust Load Performance (Bending Moment) 



/ 

— ' NOMINAL 

— — information 

. . . . . MATRIX 



Figure 33. Case 4R Short Period Damping 


1. Overall Relative Score --Coarse measure of the performance 
of the nominal controller with respect to each specification. 

_ Ideal Score - Score of Insensitive Controller 
Ideal Score - Score of Nominal Controller 


2. Normalized Performance /Range --The normalized performance 
is the average of the performance of the nominal controller 
divided by the average of the performance of the insensitive 
controller for each design specification. 


I|P!I 


_ V P 

N ^ i 
i 1C 


N ? 


^N 


i = 1, N 

N = number of eval- 
uation conditions 
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The normalized range is the range of the insensitive controller 
divided by the range of the nominal controller 


ll^li 


RSS 

(no. of criteria) 


MAX.(P^^ -MIN.(P^^) 
MAX.(P^) - MIN.(Pj^) 


i = 1, N 


3. Normalized Specification Violation- -Total spec violations for 
each insensitive controller for all evaluation conditions nor- 
malized by the maximum spec violation 


!lsv II 






RSS 

. sv. 
1 1 

i = 

1. N 

(no. of criteria) 

MAX E SV 
1. 

J = 

1, M 



M = 

no. of insensi- 
tive controllers 


The criteria are described in more detail in Reference 1. As they have been 
defined, the lower the numerical rating the better the performance of the 
controller. 


Figure 34 shows the performance of the Information Matrix controller as 
measured by the overall relative score versus the eight controllers evaluated 
in Reference 1. Figures 35 and 36 show similar comparisons for the Infor- 
mation Matrix controller performance as measured by the Normalized Per- 
formance/Range score and the Normalized Spec Violation score, respec- 
tively. Table 12 presents a summary of the rankings of the Information 
Matrix controller versus the eight controllers evaluated in Reference 1. 
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OVERALL RELATIVE SCORE 


CASE 4R 
(6 CONDITIONS) 


n CONTROL ACTIVITY ADDED 
LJ TO ABOVE CASES 



Figure 34. Overall Relative Score Comparison 
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RSS PERFORMANCE/RANGE 


D 

S 


3.5 r 


RSS PER- 
FORMANCE 

Q 

RSS RANGE 

n 

U 
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RSS PERFORMANCE 
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CONTROL ACTIVITY 
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Figure 35. Normalized Performance/ Range Comparison 



NORMALIZED SPEC VIOLATION 


3 CASE 4R (6 CONDITIONS) 
^ (WITHOUT TORSION) 



Figure 36 


Normalized Spec Violation Comparison 




TABLE 12. RANKING OF INSENSITIVE CONTROLLERS INCLUDING 
THE INFORMATION MATRIX APPROACH 


Controller 


Overall 

Relative 

Scoring 


Overall 

Performance / Range 


Overall 
Specification 
Violation S 


Information 

Matrix 

3 

Minimax 

1 

Uncertainty 

Weighting 

2 

Additive 

Noise 

4 

Multiplant 

5 

Mismatch 

Estimation 

6 

Nominal 

7 

State- 

Dependent 

Noise 

8 

Sensitivity 

Vector 

9 






Comparisons and Conclusions 


Although the Information Matrix controller does result in improved perfor- 
mance over the top-ranked minimax and uncertainty weighting of Reference 1, 
it is premature to state that the Information Matrix technique is in some way 
better than the others. The fact that the controller designed with the Infor- 
mation Matrix technique did do well, though, indicates significant potential 
for the approach. It is extremely difficult to extend a theoretical concept to 
practical worthiness with one application. It is felt that further investigation 
into increased values of would have resulted in an even better perfor- 
mance. It is also felt that a reformulation of the response vector to better 
control short period frequency and damping would have improved perfor- 
mance. (This actually refers to the evaluation of the insensitive controllers 
in Reference 1.) Some definite advantages that can be stated at this time 
include: 

1. No a priori range of parameter variations is required since 
the design is done at the nominal. One needs only the partial 
derivatives of the system matrices at the design point. 

2. The method treats nicely the response uncertainties. 

3. Since the control is assumed to be in the form u = Kx, only 
n^ -dimensional equations need be solved. No extra modeling 
or filters are necessary. 

4. It treats the actual closed-loop sensitivity, unlike the sensi- 
tivity vector augmentation or uncertainty weighting methods . 
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5. The technique can be extended easily to limited -state feedback 
and possibly observer /Kalman filter cases. 

6. With the modified Bartels -Stewart algorithm, only a modest 
eight linear equations (equivalent) need be solved per iteration 
to get the cost and gradient for p = 3. 

7. The approach provides an intuitive feel and insight to the design 

problem. It indicates clearly the improvement in system sen- 
sitivity (in terms of the dispersion matrix) and the price paid 
in terms of performance The key parameter Pg controls 

the trade-off between sensitivity and performance. 

8. The technique can be tuned to weight the relative importance 
of different parameters. 
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SECTION V 


CONCLUSIONS AND RECOMMENDATIONS 


The objective of this study, to develop useful synthesis techniques from the 
two advanced theoretical concepts created in the previous study, has been 
satisfied. This study has shown that 

• The insensitive controller synthesis technique based on the 
Finite Dimensional Inverse (FDI) concept is impractical for 
flight control system design in its current formulation. This 
is due to the time -varying nature of the resultant controller, 
which is more amenable to trajectory-type applications. 

• Despite severe computational requirements, FDI controller 
synthesis and implementation are feasible. Experiments 
with recycling stored data to alleviate storage requirements 
produced apparently satisfactory results following some 
initialization transients that could be reduced with filtering. 

• The FDI technique provides an on-line identification capability 
that could be useful for many applications. 

• The Information Matrix (IM) synthesis technique is definitely 
applicable to flight control problems because the resultant 
insensitive controller has constant gains. 

• The IM controller performs as well as the top-ranked 
uncertainty weighting and minimax controllers of the previous 
study. As in the previous study, it will be necessary to 
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qualify the results of the evaluation. It should be emphasized 
that good performance on one example with one set of criteria 
does not imply universal goodness. The IM approach, how- 
ever, with its design feature of weighting performance versus 
sensitivity without specifying the range of uncertain param- 
eters, together with the evaluation results indicate a worth- 
while development. 

Based on the results of this study, we recommend the following areas for 
further research: 

• Formulation of the FDI approach to handle stationary 
problems . 

• Development of FDI capability to handle nonminimum phase 
systems. 

This capability is needed, as demonstrated in the C-5A 
example, when the design responses used result in unstable 
compensation. 

• Further refinement of the IM methodology to quantify the 
import of modulation of controller design parameters. 
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APPENDIX A 


THE SEVENTH -ORDER MODEL 


The numerical data for the seventh-order model for six parameter values 
is given. The usual state variable representation is used: 

X = Fx + G^u + G^T] 
r = Hx + Du 

with 


T 

x 

= [w, q/n 2 , i\^. 

fia, ae^, Wg] 

T _ 

u 

= *®icJ 


T 

• • 

• • 

r 

~ ^1’ ^1* "^1* 

ta, «ej, sa^. ee., 


The Gj^ and G 2 matrices are the same for all parameter values. They are: 





r~ 

0 

0 


0 

0 

0 


0 

0 

0 


0 

0 

0 

^2 = 

0 

6 

0 


0 

0 

7.5 


0 

0 

0 


0.861 


— 


w _ 
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The last four rows of F and the last five rows of H and D are also indepen- 
dent of the parameters. Thus, we may write F, H, and D as: 



The matrices F^, and are shown below for the six parameter 

values in Tables Al through A6. 
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TABLE Al. NOMINAL CONDITION MATRICES 









Row 1 

-6.9991E-01 

3, 2724E+00 

-3.5667E-02 

-6. 6031E-01 

-2. 0358E+02 

-2. 0256E+02 

-8.6914E+00 

Row 2 

-4. 0476E-01 
Row 3 

-1, 0959E+00 

3. 1639E-02 

-1. 8709E-01 

-5. 2655E+02 

-2.3074E+03 

-9.4898E-K)0 

-1. 6852E+00 

1, 8052E-02 

-9. 7983E-01 

-3.0030E+01 

-3, 0981E+03 

1. 235 7E+03 

-1,8459E+01 








Row 1 

-2. 8733E+04 
Row 2 

-1.5687E+03 

1. 3174E+04 

1, 2412E+06 

3. 25 75E+07 

3.5344E+06 

-3.3905E+05 

-2. 2396E+04 
Row 3 

-3.4449E+03 

-1. 6286E+03 

-8.6938E+04 

8.6275E+06 

-2.4293E+06 

-2. 7210E+05 

-1.5764E+04 
j Row 4 

-1.4707E+04 

1. 2510E+06 

-3.3848E+05 

-2.5627E+07 

1, 2064E+07 

-1.4803E+05 

5.4500E+03 

-2, 7842E+04 

-5. 8879E+04 

-6.6787E+03 

1, 3327E+07 

2.6957E+06 

7, 1959E+04 








Row 1 
0. 

Row 2 

0. 






0. 

0, 






Row 3 

-8.4864E+06 
Row 4 

3. 1845E+06 






-1. 2719E+07 

2.4720E4-05 









TABLE A2. WORST CASE 1 MATRICES 


; 







Row 1 

-9.5644E-01 
1 Row 2 

4.0757E+00 

-5.9648E-02 

-1.0294E+00 

-1.6296E+02 

-2.6206E+02 

-1. 1848E+01 

-4.9662E-01 

-1. 2354E+00 

4.0758E-03 

-5.9496E-01 

-5, 2179E+02 

-2.5402E+03 

-1.1194E+01 

Row 3 

-2.9321E+00 

-2. 1624E-01 

-1, 2407E+00 

-2.394E+01 

-2.9521E+03 

1.2809E+03 

-3.3386E+01 








Row 1 

-3. 7187E+04 

-3, 1421E+03 

1.5809E+04 

1.4645E+06 

3.9282E+07 

1.5362E+06 

-4.4309E+05 

Row 2 

-3.0729E+04 
Row 3 

-5.0614E+03 

-2. 1960E+03 

-7. 0371E+04 

1.3523E+07 

-4.0092E+06 

-3. 7416E+05 

-3.4336E+04 
Row 4 

-2.6231E+04 

1.4743E+06 

-3, 2140E+05 

-2. 8624E+07 

1.6160E+07 

-3.6035E+05 

8. 1625E+03 

-4.3576E+04 

-2.9472E+04 

-7.6219E+03 

1.5514E+07 

4.4698E+06 

1.0772E+05 








Row 1 

0. 

Row 2 

0. 






0. 

0. 






Row 3 

: -1.0608E+07 

Row 4 

3.9806E+06 






-1.5899E+07 

3.0900E+05 
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TABLE A5, uj PERTURBATION MATRICES 
(tt) = 0.75 









Row 1 

-7.3843E-01 
Row 2 

3. 2656E+00 

-4. 2888E-02 

-7.5762E-01 

-1, 6024E+02 

-2. 0636E+02 

-9. 1547E+00 

-5. 0739E-01 
Row 3 

-1. 0219E+00 

1. 3385E-02 

-3.6874E-01 

-4.5964E+02 

-2. 1215E+03 

-9.0558E+00 

-2, 0814E+00 

-1.0020E-01 

-9. 8431E-01 

-2. 1095E+01 

-2.6579E+03 

1, 1008E+03 

-2.3423E+01 








Row 1 

-2.9371E+04 

-2. 1935E+03 

1. 2852E+04 

1. 1860E+06 

3. 179 7E+07 

2. 005 2E+06 

-3.4887E+05 

Row 2 

-2.3701E+04 
Row 3 

-3. 8144E+03 

-1. 6975E+03 

-6. 1585E+04 

9.9185E+06 

-2. 9214E+06 

-2. 8838E+05 

-1.9502E+04 
Row 4 

-1.6019E+04 

1. 1962E+06 

-2. 3069E+05 

-2.0137E+07 

1.0730E+07 

-1.9624E+05 

5.3642E+03 

-2. 7868E+04 

-3. 0736E+04 

-4. 9718E+03 

1. 2773E+07 

2. 76 75E+06 

7.0568E+04 








Row 1 

0. 

Row 2 

0. 






0. 

0. 






Row 3 

-8.4864E+06 
Row 4 

3. 1845E+06 






-1. 2719E+07 

2. 4720E+05 
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TABLE A6, q PERTURBATION MATRICES 
(q = 1- 25q^Q^) 









Row 1 

-8.8911E-01 
Row 2 

4.0881E+00 

-4. 7302E-02 

-8.6178E-01 

-2.3841E+02 

-2.5439E+02 

-1. 1035E+01 

-6.4055E-01 
Row 3 

-1.3372E+00 

3. 2041E-02 

-3.0658E-01 

-6.3226E+02 

-2.8037E+03 

-1. 1641E+01 

-2.2559E+00 

-2. 3891E-02 

-1. 1964E+00 

-3, 2267E+01 

-3. 7075E+03 

1.4896E+03 

-2.4954E+01 








Row 1 

-3.6171E+04 
Row 2 

-2, 2256E+03 

1.6183E+04 

1.5200E+06 

4.0385E+07 

3. 7692E+06 

-4. 2780E+05 

-2.8484E+04 
Row 3 

-4.4490E+03 

-1.9753E+03 

-9,3611E+04 

1. 1260E+07 

-3.2345E+06 

-3.4624E+05 

-2.6015E+04 
Row 4 

-2.3634E+04 

1.5296E+06 

-4.4507E+05 

-3.9464E+07 

1.9481E+07 

-2.5615E+05 

8.3336E+03 

-4.3517E+04 

-5. 7628E+04 

-1. 1698E+04 

1.6277E+07 

4.4470E+06 

1.0918E+05 








Row 1 

0. 

Row 2 

0. 






0. 

0. 






Row 3 

-1.0608E+07 
Row 4 

3.9806E+06 






-1.5899E+07 

3.0900E+05 








APPENDIX B 


MAXIMUM DIFFICULTY METRIC 


This appendix summarizes work performed (on the maximum difficulty 
metric) in determining the design point in parameter space at which control 
is most difficult. 

Recall that the basic idea was to find the point a in parameter space to 
maximize the difficulty metric 

J = trlH' (H W HO'^ H] = tr D (B-1) 

a 

The approach taken is to get an analytic expression for the gradient 5>J/5a 
that could be used as a basis for a numerical optimization scheme. 

In Equation (B-1) 

A AT 

H = He (B-2a) 

T -Af -1 -A 

W = e BR B"e dt (B-2b) 

where T is an arbitrary parameter. Also, for notational convenience 

^ A ^ ^ 

X = (H W H') 

The matrices A, H are subject to parameter uncertainty 
p 

A-.A + Z A. (B-3a) 

i=l 
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(B-3b) 



P 

H -♦ H + 2 

Consider now a first order pertubation in -» and the resulting 

change in J: 

J + 6J = tr[(H + 6 H) ' {X + 6 X}"^ (H + 6 H)] 
expanding to first order, keeping terms of o( 6) only and using 

(X+ ^ X"^ - X"^ (6X) X"^ 

we obtain 

«J = 2 tr (H'X“^ 5H) - tr (X"^ H H"X fiX) (B-4) 

A 

So we need only to get flH and sX in terms of 


First consider 6H. 


A AT* At* 

6H = (5H) e^^ + H (6e^^) 


XT AT (A + Scx.A.)T AT ^ 

= 6a. H. e + H fe i i -e 1 

1 1 


Using a result from Bellman, 
(A + Sof. A.)T AT 


1 1 


, ^ „T A(T - t) . AT 
+ 6a. e A. e dt 

1 >J o 1 


AT ^ AT 
= e + 6aj^ e r (T) 


Thus, 


^ = IS = H. e^"^ + H r (T) 


(B-5) 

(B-5a) 

(B-6) 
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where 


(»T -AT . AT 
r(T) = e A^ e dt 


(B-7) 


The process of computing 8X proceeds in a similar manner. 


W + 6W = r'^ BR-l B' ‘“i ^i>'‘ 


^ o 


dt 


It is necessary to work with X, so with a change of variable, 

X + 6X = (H + 60 '. H.) [r"^ ®"i "^i^^ BR"^B'e^^ ^ 6a. A.) 'cj^ ^ 

— — 11 '-'Jo ■’ 

(H -1- 6a. H.)' 


Expanding the integral term, using Equation (B-5), gives 

'T A _ _1 A /. 


— = H [(’’^ e^"^ BR"^ B' r' (cj) da + f ^ e^'^ r(a)BR"-'B^e^ ^da]H' 

tot. O " o 


T Ac 


-1 


A A ^T AT 

+ H W e H.' -t- H.e W H 
1 1 


Substituting Equations (B-6) and (B-8) into Equation (B-4) gives 

^ = 2 tr (H X — ) - tr (X H H"X ~) 

5of. - ba. - - ba^ 

1 At 2 1 AT 

= 2tr[H'X H.e 4-Dr(T)-D W.,(T)-D HX H. e W1 

— 1 Of orl O'— 1 


* H W H " = H f’ ■" e 
^ O 


T „A(T-t) '^'^■^^dtH ' = H l’'^ e'^'^BR"^B' e'^^'^daH ' 


- o 
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where D is the difficulty matrix and 


a 


W, (T) = e'^^ BR”^ B' r' (T - t) dt 

1 - O 


J’j' BR“’- B' e“^'^ [J^ A.^ e'^'^ d 5 ] dt 


Computational methods for evaluating Wj^(T) have been considered, but 
thus far nothing simple has come up. A straightforward numerical evaluation 

_ J (or + 5^.) - J (a ~ Sof.) 
i 


may be easier. 
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